Hyperparameter Search and Nested Cross-Validation › Formula

It helps to view hyperparameter search as a second learning problem wrapped around the first. The inner problem fits model parameters for a fixed hyperparameter setting. The outer problem compares those settings using held-out data. The formulas below make both layers explicit, because the distinction between them is exactly what nested cross-validation is designed to protect.

Let the training data be

$$D=\{(x_i,y_i)\}_{i=1}^{n},$$

where $x_i$ is a feature vector and $y_i$ is a target. A learning algorithm $A$ takes a training set and a hyperparameter vector $\lambda$ and returns a fitted model $f_\lambda$. The hyperparameters live in a search space $\Lambda$, where:

• $\lambda$: a hyperparameter vector, an element of the search space
• $\Lambda$: the search space, the set of all allowed hyperparameter values
• $A$: the learning algorithm that maps a training set and $\lambda$ to a fitted model
• $f_\lambda$: the fitted model returned by $A$ for hyperparameters $\lambda$

In a polynomial ridge pipeline, for example,

$$\lambda=(d,\alpha),$$

where:

• $d$: the polynomial degree (a discrete hyperparameter)
• $\alpha$: the ridge penalty strength (a continuous hyperparameter)

The degree controls the feature map

$$\phi_d(x)=(x,x^2,\ldots,x^d),$$

where:

• $\phi_d(x)$: the polynomial feature map that sends a scalar input $x$ to a vector of its powers up to degree $d$
• $d$: the polynomial degree, the highest power included in the feature map
• $x$: a scalar input feature value

Note that $\phi_d(x)$ is a $d$-dimensional vector for scalar input $x$, with components $x, x^2, \ldots, x^d$, so the coefficient vector $\beta$ also has length $d$, ensuring the inner product is well defined:

$$\phi_d(x_i)^\top\beta = \sum_{p=1}^{d} \beta_{p} x_i^{p}$$

The dot product is a scalar prediction, where:

• $\phi_d(x_i)$: the $d$-dimensional polynomial feature vector for observation $i$, with components $x_i, x_i^2, \ldots, x_i^d$
• $\beta$: the coefficient vector of length $d$ that weights each power
• $\beta_p$: the coefficient associated with the $p$-th power of the input
• $x_i$: the feature value for observation $i$ (a scalar in this univariate polynomial case)
• $d$: the polynomial degree, the highest power included in the sum

This feature map excludes a constant column. The implementation adds an unpenalized intercept and standardizes each power using training-fold means $\mu_p$ and positive scales $s_p$. Write the standardized power as $z_{ip}=(x_i^p-\mu_p)/s_p$ and the standardized feature vector as $z_i=(z_{i1},\ldots,z_{id})$, where:

• $z_{ip}$: the standardized version of the $p$-th power of input $x_i$
• $\mu_p$: the mean of $x^p$ computed on the training fold only
• $s_p$: the positive scale of $x^p$ (typically its standard deviation) computed on the training fold only
• $z_i$: the vector of standardized polynomial features for observation $i$

For a constant column, use scale one. Ridge trades residual error against coefficient size on this standardized scale:

$$\left(\hat\beta_0(d,\alpha),\hat\beta_\alpha\right)=\arg\min_{\beta_0,\beta}\left(\sum_{i=1}^{n}\left(y_i-\beta_0-z_i^\top\beta\right)^2+\alpha\lVert\beta\rVert_2^2\right),$$

where:

• $\hat\beta_0(d,\alpha)$: the unpenalized intercept estimate, a function of the degree $d$ and penalty $\alpha$
• $\hat\beta_\alpha$: the fitted (penalized) coefficient vector at penalty strength $\alpha$; it implicitly depends on the standardized feature basis determined by $d$
• $\beta$: the coefficient vector being optimized over, of dimension $d$ (the same length as the standardized feature vector $z_i$)
• $n$: the number of training observations
• $y_i$: the observed target for observation $i$
• $\beta_0$: the unpenalized intercept, not penalized by the ridge term
• $z_i$: the polynomial features standardized using only the training observations
• $\beta_0+z_i^\top\beta$: the model's predicted value for observation $i$
• $\alpha$: the ridge penalty strength (controls the size of the shrinkage)
• $\lVert\beta\rVert_2^2$: the squared L2 norm of the coefficients, $\sum_{p=1}^{d}\beta_{p}^2$

We will study linear regression and polynomial regression in detail later; for now, the important point is that $d$ changes the model family and $\alpha$ changes the fitting objective, but neither is learned by the inner least-squares step. Changing $\lambda$ changes the model family or the fitting rule. Note that the ridge objective above shows the unpenalized intercept $\beta_0$ alongside the penalized coefficients, so the prediction used in the loss includes both.

The quantity we would like to know is the generalization risk

$$R(\lambda\mid D)=\mathbb{E}_{(X,Y)\sim P}\left[L\left(Y,\,A(D;\lambda)(X)\right)\right],$$

where:

• $R(\lambda\mid D)$: the expected loss of the particular model fitted on $D$ with setting $\lambda$, where $\lambda$ is a hyperparameter vector in the search space $\Lambda$
• $D$: the training data, held fixed in this conditional risk
• $(X,Y)$: a new observation drawn from the data-generating distribution
• $P$: the data-generating distribution
• $L$: a loss such as squared error
• $A(D;\lambda)$: the fitted model returned by the learning algorithm on training data $D$ with setting $\lambda$
• $X$: a new input feature vector drawn from the data-generating distribution
• $A(D;\lambda)(X)$: the fitted model's prediction at the new input $X$

and the squared-error loss is

$$L(y,\hat y)=(y-\hat y)^2,$$

where:

• $L(y,\hat y)$: the squared-error loss for a true target $y$ and prediction $\hat y$
• $y$: the observed target
• $\hat y$: the model's predicted value

The new observation must be independent of the training data. In this conditional expression, fitting randomness is fixed; the expectation averages over new evaluation observations. A procedure-level risk additionally averages this conditional risk over fresh training samples of a specified size and over any fitting randomness. These are different targets: cross-validation evaluates fits trained on fewer rows than a fit on all $n$ observations. Below, $R(\lambda)$ abbreviates the procedure-level risk at the training size being evaluated, with that size held fixed across candidates. We estimate risk with held-out data, once per candidate.

For $K$-fold cross-validation, partition the training indices into folds $V_1,\ldots,V_K$. Let $D^{(-k)}$ be the data with fold $V_k$ removed, and let

$$f_\lambda^{(-k)}=A\left(D^{(-k)};\lambda\right)$$

be the model fitted without fold $k$. The cross-validation estimate averages the loss over all folds and, within each fold, over the held-out observations, so that every observation contributes once to the estimate:

$$\hat R_{\mathrm{CV}}(\lambda) = \frac{1}{K}\sum_{k=1}^{K} \frac{1}{|V_k|} \sum_{i\in V_k} L\left(y_i,f_\lambda^{(-k)}(x_i)\right)$$

where:

• $\hat R_{\mathrm{CV}}(\lambda)$: the $K$-fold cross-validated estimate of the risk at hyperparameters $\lambda$
• $K$: the number of folds
• $k$: the index of the fold currently held out
• $V_k$: the index set of observations in fold $k$
• $|V_k|$: the number of observations in fold $k$
• $f_\lambda^{(-k)}$: the model fitted on all folds except $k$, with hyperparameters $\lambda$
• $L$: the per-observation loss
• $y_i$: the observed target for observation $i$
• $x_i$: the feature vector for observation $i$

The inner sum averages loss within held-out fold $k$, and the outer sum gives each fold equal weight. If fold sizes differ, this need not equal the pooled per-observation mean loss. The superscript $(-k)$ is important: the model used to predict observations in $V_k$ never saw those observations during fitting. For a fixed $\lambda$, this is a standard validation estimate. The subtlety arises because we will compare many such estimates and keep the smallest one.

Parameters are the numbers that the fitting algorithm optimizes. In a linear model, they are the regression coefficients; in a neural network, they are the weights and biases. Hyperparameters are the numbers or structural choices that control how that optimization is set up. They include the polynomial degree, the ridge penalty, the number of folds in cross-validation, the learning rate of gradient descent, the number of trees in a forest, and the kernel width in a kernel method. The distinction is practical: parameters are estimated by minimizing a training objective, while hyperparameters are chosen by comparing validation performance or by using prior knowledge.

A useful way to keep the distinction straight is to ask what changes when a value changes:

• A parameter changes during fitting on a fixed training set. The fitting rule itself stays the same.
• A hyperparameter changes the fitting rule, the model family, or the search protocol. It is fixed while the parameters are optimized.
• A search space specifies which hyperparameter values are allowed. It is part of the procedure and must be fixed before the final evaluation.

The search space $\Lambda$ can contain discrete values, continuous intervals, categorical choices, or conditional choices. For a polynomial ridge model we might use the degree and the penalty:

$$d\in\{1,2,\ldots,8\}, \qquad \alpha\in[10^{-3},1],$$

where:

• $d$: the polynomial degree (a discrete hyperparameter)
• $\alpha$: the ridge penalty strength (a continuous hyperparameter)

For a positive penalty spanning several orders of magnitude, a logarithmic search allocates resolution to relative changes in $\alpha$. Moving from $0.001$ to $0.01$ changes the penalty by a factor of ten, just as moving from $0.1$ to $1$ does. The interval $[0,0.01]$ occupies only 1% of the linear range $[0,1]$. A coarse linear grid can therefore miss substantial relative differences among small positive penalties. Search spaces should be chosen before looking at the final test result, and they should be wide enough to include reasonable models but not so wide that the budget is spread too thin.

Once the search space is defined, the next question is how to spend the evaluation budget. Grid search evaluates a Cartesian product of candidate values. If the search space has $H$ hyperparameters and the $j$-th hyperparameter has $m_j$ candidate values, the grid is

$$\Lambda_{\mathrm{grid}} = \left\{ \left(\lambda_1^{(a_1)},\ldots,\lambda_H^{(a_H)}\right) : a_1\in\{1,\ldots,m_1\},\ldots,a_H\in\{1,\ldots,m_H\} \right\},$$

where:

• $\lambda_j^{(a_j)}$: the $a_j$-th candidate value for the $j$-th hyperparameter
• $m_j$: the number of candidate values for the $j$-th hyperparameter
• $H$: the number of hyperparameters being tuned
• $a_j$: the index selecting the candidate value along the $j$-th hyperparameter dimension

Each element of $\Lambda_{\mathrm{grid}}$ is a single combination, chosen by fixing one candidate index $a_j$ per dimension. The number of candidates is the size of the Cartesian product, obtained by multiplying the number of candidate values along each hyperparameter dimension:

$$B=\prod_{j=1}^{H} m_j,$$

where:

• $B$: the total number of grid candidates
• $H$: the number of hyperparameters being tuned
• $m_j$: the number of candidate values for the $j$-th hyperparameter
• $\prod_{j=1}^{H}$: the product over all $H$ hyperparameters
• $j$: the index identifying the hyperparameter whose candidate count is multiplied

If we tune eight polynomial degrees and four ridge penalties, the budget is $8\times 4=32$ candidate configurations. Each configuration is evaluated on the same inner folds. Grid candidate enumeration is deterministic and exhaustive over the specified grid. Randomness in model fitting or fold construction still needs to be controlled. Its weakness is the product: five hyperparameters with ten values each require $10^5=100{,}000$ candidates, which is often impractical. The grid also evaluates every combination, including combinations whose differences may be irrelevant to performance.

Random search draws candidates independently from a distribution $p(\lambda)$:

$$\lambda^{(b)}\sim p(\lambda), \qquad b=1,\ldots,B.$$

where:

• $\lambda^{(b)}$: the $b$-th candidate hyperparameter vector drawn by the search
• $p(\lambda)$: the search distribution over the hyperparameter space
• $B$: the total number of candidates (the budget)
• $b$: the index of the candidate draw

We then select

$$b^*=\arg\min_{b\in\{1,\ldots,B\}}\hat R_{\mathrm{CV}}\left(\lambda^{(b)}\right), \qquad \hat\lambda=\lambda^{(b^*)}.$$

where:

• $b^*$: the selected candidate index, using the first index to break an exact tie
• $\hat\lambda$: the selected candidate vector at that index
• $\hat R_{\mathrm{CV}}\left(\lambda^{(b)}\right)$: the $K$-fold cross-validation estimate of the risk for the $b$-th candidate

For the positive ridge penalty in this chapter, we choose a log-uniform distribution because it gives equal probability to equal-width intervals on the logarithmic scale, including each full order of magnitude within its support. Random search does not enumerate combinations. Independent draws often produce many different values along continuous dimensions; discrete values can repeat, and distinctness is not guaranteed. If only one or two hyperparameters strongly affect performance, random search can cover their important ranges much faster than a grid that also varies irrelevant dimensions. It is stochastic, however, so the result depends on the seed, and two runs with the same budget can select different configurations.

An equal-budget comparison must hold the data, inner folds, loss, and number of candidate evaluations fixed. Grid search and random search should not be compared when one method receives more fits or different validation splits. For random search, we should report the seed or repeat the search across several seeds, because the candidate draw is part of the procedure. The fair question is not "which method found a lower score once?" but "how does each complete search procedure behave under the same budget and protocol?"

Let $\Lambda_B$ denote the finite set of evaluated candidates within the larger allowed space $\Lambda$. After evaluating those candidates, we choose the one with the smallest validation loss:

$$\hat\lambda = \arg\min_{\lambda\in\Lambda_B} \hat R_{\mathrm{CV}}(\lambda).$$

The reported score of the winner is therefore

$$\hat R_{\mathrm{CV}}(\hat\lambda) = \min_{\lambda\in\Lambda_B} \hat R_{\mathrm{CV}}(\lambda).$$

This is a minimum of noisy estimates. Even if all candidates have the same true risk, the one that happens to receive a favorable validation sample will look better than the others. To see the mechanism, write the validation error of candidate $\lambda$ as

$$\hat R_{\mathrm{CV}}(\lambda) = R(\lambda)+\varepsilon_\lambda,$$

where:

• $\hat R_{\mathrm{CV}}(\lambda)$: the cross-validated estimate of the risk at hyperparameters $\lambda$
• $R(\lambda)$: the true risk of the model fitted with hyperparameters $\lambda$
• $\varepsilon_\lambda$: the estimation error on these particular folds (the gap between the validation estimate and the true risk)

The selected candidate has

$$\hat R_{\mathrm{CV}}(\hat\lambda) = R(\hat\lambda)+\varepsilon_{\hat\lambda},$$

where:

• $\hat R_{\mathrm{CV}}(\hat\lambda)$: the validation estimate of the selected candidate's risk
• $R(\hat\lambda)$: the true risk of the selected candidate
• $\varepsilon_{\hat\lambda}$: the estimation error of the selected candidate on these folds, the gap between its validation estimate and its true risk (defined as $\hat R_{\mathrm{CV}}(\lambda)-R(\lambda)$ evaluated at $\hat\lambda$)

Selection favors candidates with favorable estimation errors. Under the unbiased-candidate assumptions stated below, the resulting optimism is nonnegative in expectation; this is not a claim about every realized score. Defining optimism as the expected amount by which the selected model's validation estimate falls below its true risk (equivalently, the expected overstatement of performance), we have

$$\begin{aligned} \mathrm{optimism} &= \mathbb{E}\left[R(\hat\lambda)\right] - \mathbb{E}\left[\hat R_{\mathrm{CV}}(\hat\lambda)\right] \\ &= -\mathbb{E}\left[\varepsilon_{\hat\lambda}\right] \end{aligned}$$

where:

• $\mathrm{optimism}$: the expected overstatement of performance caused by selecting the best candidate
• $\hat\lambda$: the selected hyperparameters (a random quantity, depending on the validation folds)
• $R(\hat\lambda)$: the true risk of the model chosen by the search
• $\hat R_{\mathrm{CV}}(\hat\lambda)$: the validation estimate of the selected model's risk
• $\mathbb{E}[\cdot]$: the expectation over the randomness of the validation folds

For a fixed finite candidate set, suppose each candidate estimate is unbiased for its fixed target risk: $\mathbb{E}[\hat R_{\mathrm{CV}}(\lambda)]=R(\lambda)$. The minimum is a concave function of its arguments. Applying Jensen's inequality to this concave function gives $\mathbb{E}[\min_\lambda \hat R_{\mathrm{CV}}(\lambda)]\le\min_\lambda \mathbb{E}[\hat R_{\mathrm{CV}}(\lambda)]$, and substituting the unbiasedness $\mathbb{E}[\hat R_{\mathrm{CV}}(\lambda)]=R(\lambda)$ yields $\mathbb{E}[\min_\lambda\hat R_{\mathrm{CV}}(\lambda)]\le\min_\lambda R(\lambda)$. Also, over the fixed finite candidate set $\Lambda_B$, $\min_{\lambda\in\Lambda_B} R(\lambda)\le R(\hat\lambda)$ for every realization (since $\hat\lambda\in\Lambda_B$), so $\min_{\lambda\in\Lambda_B} R(\lambda)\le\mathbb{E}[R(\hat\lambda)]$. Combining these gives $\mathbb{E}[\min_{\lambda\in\Lambda_B}\hat R_{\mathrm{CV}}(\lambda)]\le\mathbb{E}[R(\hat\lambda)]$, so the expected optimism is nonnegative. Independence across candidates is unnecessary for this inequality. Actual CV estimates need not be unbiased for the risk of a model refitted on all data; the training-size difference is another source of discrepancy. The inequality is not a guarantee that a realized outer score exceeds its inner score.

A simple calculation makes the size of the effect concrete. Suppose two candidates have the same true risk, and their validation estimates are independent normal random variables with the same mean $R_0$ and positive standard deviation $\sigma$. The normality and independence assumptions are stated explicitly because the exact value below depends on them; normality is an illustrative model assumption rather than a general property of validation estimates, and shared folds can make candidate estimates correlated. Write the standardized estimates as $Z_1, Z_2 \sim \mathcal{N}(0,1)$, independently. The expected minimum is

$$\begin{aligned} \mathbb{E}[\min(Z_1, Z_2)] &= \int_{-\infty}^{\infty} z \cdot 2\,\varphi(z)\,\Phi(-z)\, dz \\ &= -\frac{1}{\sqrt{\pi}} \\ &\approx -0.564 \end{aligned}$$

where:

• $z$: the value at which we evaluate whether one draw is the smaller
• $\varphi(z)$: the standard normal density
• $\Phi(z)$: the standard normal cumulative distribution function, the probability that a standard normal draw is at most $z$; its exceedance probability is $1-\Phi(z)=\Phi(-z)$
• the factor $2\,\varphi(z)\,\Phi(-z)$: the density that $z$ is the smaller of the two draws (one draw at $z$, the other above it, in either order)

To evaluate the integral, use the identity $\mathbb{E}[\min(Z_1, Z_2)] = -\mathbb{E}\lvert Z_1 - Z_2\rvert / 2$. With $W = Z_1 - Z_2 \sim \mathcal{N}(0, 2)$, we have $\sigma_W = \sqrt{2}$ and $\mathbb{E}\lvert W\rvert = \sqrt{2/\pi}\,\sigma_W$ for $W \sim \mathcal{N}(0, \sigma_W^2)$. This follows from the half-normal expectation

$$\begin{aligned} \mathbb{E}|W| &= 2\int_0^\infty w \cdot \frac{1}{\sqrt{2\pi}\sigma_W} e^{-w^2/(2\sigma_W^2)} \, dw \\ &= \sqrt{\frac{2}{\pi}}\,\sigma_W && \text{(half-normal mean for } W\sim\mathcal{N}(0,\sigma_W^2)\text{)} \\ &= \sqrt{\frac{2}{\pi}}\,\sqrt{2} && \text{(substitute } \sigma_W=\sqrt{2}\text{)} \\ &= \frac{2}{\sqrt{\pi}}. \end{aligned}$$

where:

• $w$: the integration variable running over nonnegative values
• $\sigma_W$: the standard deviation of the zero-mean normal $W$ (equal to $\sqrt{2}$ here, since $W = Z_1 - Z_2$)
• $e^{-w^2/(2\sigma_W^2)}$: the Gaussian kernel evaluated at $w$, with the normalizing factor $1/(\sqrt{2\pi}\sigma_W)$ outside
• the factor $2$: accounts for the two symmetric halves of the normal density when integrating only over $w \ge 0$

Therefore

$$\begin{aligned} \mathbb{E}[\min(Z_1, Z_2)] &= -\frac{\mathbb{E}|Z_1 - Z_2|}{2} \\ &= -\frac{1}{2} \cdot \sqrt{\frac{2}{\pi}} \cdot \sqrt{2} \\ &= -\frac{1}{2} \cdot \frac{2}{\sqrt{\pi}} \\ &= -\frac{1}{\sqrt{\pi}}. \end{aligned}$$

Rescaling back to the candidates' units, the winner's validation score is about $0.564\sigma$ below the common mean even though neither candidate is better. With more candidates, the minimum is pushed even lower. More precisely, if $M$ candidates share the same true risk and their validation estimates have standard deviation $\sigma$, the expected minimum score is

$$R_0+\sigma\, \mathbb{E}[\min(Z_1,\ldots,Z_M)],$$

where:

• $R_0$: the common true risk shared by all $M$ candidates
• $M$: the number of candidates sharing that common true risk
• $\sigma$: the standard deviation of each candidate's validation estimate
• $Z_m$: the standardized validation error of candidate $m$, defined as $Z_m=(\hat R_m-R_0)/\sigma$, an independent standard normal draw

This expression is obtained by writing each estimate as $\hat R_m = R_0 + \sigma Z_m$ with $Z_m=(\hat R_m-R_0)/\sigma$ standard normal, so the minimum is $R_0 + \sigma \min_m Z_m$ and linearity of expectation gives the expression above. Since $\mathbb{E}[\min(Z_1,\ldots,Z_M)]<0$ for $M\ge 2$ under these assumptions, the expected minimum lies below the common risk $R_0$.

This expression decreases as $M$ grows because the standardized minimum $\mathbb{E}[\min(Z_1,\ldots,Z_M)]$ decreases with $M$. This is the same "winner's curse" that can arise when we select the maximum of noisy measurements.

The winner's curse can be quantified without fitting a model. For $M$ independent candidate estimates with the same true risk, the expected best observed score is $R_0$ plus $\sigma$ times the expected minimum of $M$ standard normal draws. The next cell approximates that expected minimum by numerical integration over $[-8,8]$, and the figure plots the negative value as expected selection optimism in standard-deviation units.

Nested cross-validation addresses this by adding an outer loop. Partition the data into outer folds $V_1^{(o)},\ldots,V_J^{(o)}$, and let $D^{(-j)}$ denote the data with outer fold $V_j^{(o)}$ removed. For each outer fold $j$, run the entire inner search on $D^{(-j)}$ and obtain

$$\hat\lambda_j = \arg\min_{\lambda\in\Lambda_B} \hat R_{\mathrm{CV}}^{(j)}(\lambda),$$

where:

• $\hat\lambda_j$: the hyperparameters selected by the inner search when outer fold $j$ is held out
• $\hat R_{\mathrm{CV}}^{(j)}(\lambda)$: the $K$-fold inner cross-validation estimate of candidate $\lambda$ computed only on $D^{(-j)}$

For random search, $\Lambda_B$ here means the candidates drawn for this outer training run. Candidate generation must also exclude the outer test targets.

Then fit the selected model on the outer training set,

$$f_j=A\left(D^{(-j)};\hat\lambda_j\right),$$

and evaluate it on the outer test fold. The nested estimate averages the outer-fold losses, each computed with the model whose hyperparameters were selected by that fold's inner search:

$$\hat R_{\mathrm{nested}} = \frac{1}{J}\sum_{j=1}^{J} \frac{1}{|V_j^{(o)}|} \sum_{i\in V_j^{(o)}} L\left(y_i,f_j(x_i)\right)$$

where:

• $\hat R_{\mathrm{nested}}$: the nested cross-validation estimate of the procedure's risk
• $J$: the number of outer folds
• $D^{(-j)}$: the data with outer fold $j$ removed, used to run the inner search and select $\hat\lambda_j$
• $V_j^{(o)}$: the index set of observations in outer fold $j$
• $|V_j^{(o)}|$: the number of observations in outer fold $j$
• $f_j$: the model refitted on the outer training set using the selected hyperparameters $\hat\lambda_j$
• $L$: the per-observation loss
• $y_i$: the observed target for observation $i$
• $x_i$: the feature vector for observation $i$
• $i$: the index running over observations within an outer test fold
• $j$: the index of the outer fold

The outer test fold was not used to choose $\hat\lambda_j$, so it gives a cleaner view of how the selection procedure performs on unseen data. Notice that $\hat\lambda_j$ can differ from fold to fold. The nested estimator is not the cross-validation score of one fixed model; it is an estimate of the risk of the complete procedure that takes a training set, runs a search, and returns a fitted model. With independent identically distributed observations and a procedure fixed before inspecting outer results, outer held-out losses evaluate that procedure at the outer-training sample size, approximately $(1-1/J)n$. They do not directly estimate the conditional risk of the final model refitted on all $n$ rows. Finite data and overlapping outer training sets still create variability. To see why the distinction matters, note that $\hat R_{\mathrm{nested}}$ averages over outer folds each of which selected its own $\hat\lambda_j$; keeping each outer test fold out of its entire selection process is what avoids reusing its scores for selection, but it also means the nested score is not the score of any single deployable model.

After nested cross-validation, we usually want a final model for deployment. The standard recipe is to run the search on all available development data, obtain $\hat\lambda$, and refit

$$\hat f=A(D;\hat\lambda).$$

The nested estimate describes the procedure; the final model is the procedure applied to all development data. If an untouched test set exists, we evaluate $\hat f$ on it once. That test score is often the cleanest estimate of final performance, and in that situation nested cross-validation is not mandatory. Nested cross-validation is the appropriate tool when we need to estimate the procedure using the same data that were available for development, or when we want to compare several complete procedures under a common protocol.

Reporting should include the search method, the search space, the budget, the inner and outer fold counts, the metric, the random seeds, and the selected hyperparameters across outer folds. A single number hides too much. It is useful to report the mean outer score, the standard deviation across outer folds, and a standard error across repeated nested runs. If we repeat the entire nested procedure $M\ge2$ times with different outer splits, we can summarize the repetition-level scores with a mean and an estimated standard error. The mean of the per-repetition nested estimates is

$$\bar R=\frac{1}{M}\sum_{m=1}^{M}\hat R_m,$$

and the estimated standard error of that mean is

$$\widehat{\mathrm{SE}}(\bar R) = \frac{1}{\sqrt{M}} \sqrt{\frac{1}{M-1} \sum_{m=1}^{M} \left(\hat R_m-\bar R\right)^2}$$

where:

• $M$: the number of repeated nested runs (each with a fresh outer split)
• $m$: the index of the repetition
• $\hat R_m$: the nested estimate produced by the $m$-th repetition
• $\bar R$: the mean nested estimate across repetitions
• $\widehat{\mathrm{SE}}(\bar R)$: the estimated standard error of $\bar R$ across repetitions
• $\left(\hat R_m-\bar R\right)^2$: the squared deviation of repetition $m$'s estimate from the mean
• $\sqrt{\frac{1}{M-1}\sum_{m=1}^{M}\left(\hat R_m-\bar R\right)^2}$: the sample standard deviation of the repetition-level scores

The formula is the sample standard deviation of the repetition-level scores divided by $\sqrt{M}$: the inner square root is the sample standard deviation across repetitions, and dividing by $\sqrt{M}$ converts that spread into a standard error for their mean.

With independently drawn partition and search seeds, this is a Monte Carlo standard error for the mean over that randomization, conditional on this fixed dataset. It can decrease as we add repetitions without adding observations. The repetition-level standard deviation describes split/search variability; neither quantity measures how results vary over new datasets or supplies a population-risk confidence interval. If the selected degree or penalty changes drastically from fold to fold, the procedure is unstable, and a simpler model or more data may be preferable.

Grid and random search treat validation loss as a black-box objective and do not use derivatives. A discrete degree cannot be updated by an ordinary continuous gradient. Some continuous hyperparameters, including ridge penalties, admit differentiable optimization under suitable conditions, but that is a different method from the searches implemented here. A validation surface need not be convex, and finite-data noise can affect which candidate wins.

Random search has a simple coverage property. The following formula makes the probability of covering the good region explicit for an independent-draw search. Suppose the good region has probability mass $q$ under the chosen search distribution. Each independent draw then has probability $q$ of producing an acceptable candidate; the value of $q$ is simply the probability mass the search distribution assigns to the good region, and a uniform distribution over the search region makes it equal to the good region's volume (or length) divided by the total volume (or length). The probability that $B$ independent draws contain at least one good candidate is

$$P(\text{at least one hit})=1-(1-q)^B,$$

where:

• $P(\text{at least one hit})$: the probability that at least one of the $B$ candidates falls in the good region
• $q$: the per-draw probability of an acceptable candidate under the search distribution
• $B$: the number of independent candidate draws
• $(1-q)^B$: the probability that all $B$ independent draws miss the good region, since each draw misses with probability $1-q$ and draws are independent

For a target hit probability $p$ with $0<p<1$ and per-draw hit probability $q$ with $0<q<1$, the smallest integer budget giving hit probability at least $p$ is

$$B_{\min}=\left\lceil\frac{\log(1-p)}{\log(1-q)}\right\rceil,$$

where:

• $B_{\min}$: the smallest integer number of independent draws needed to reach hit probability at least $p$, obtained by rounding the real-valued lower bound up to the next integer
• $p$: the target lower bound on the probability of at least one hit (with $0<p<1$)
• $q$: the probability mass of the good region under the search distribution (with $0<q<1$)
• $\lceil\cdot\rceil$: the ceiling function, which rounds its argument up to the smallest integer greater than or equal to it

To derive the budget bound, start from the requirement that the hit probability reach $p$ and rearrange the steps one at a time:

$$\begin{aligned} 1-(1-q)^B &\ge p && \text{(require hit probability at least } p \text{)} \\ (1-q)^B &\le 1-p && \text{(rearrange)} \\ B\log(1-q) &\le \log(1-p) && \text{(take natural log of both sides; } \log \text{ is increasing)} \\ B &\ge \dfrac{\log(1-p)}{\log(1-q)} && \text{(divide by } \log(1-q)<0 \text{, reversing the inequality)} \end{aligned}$$

Because $1-q\in(0,1)$, we have $\log(1-q)<0$, so the division used in the last step reverses the inequality direction. Because $1-p\in(0,1)$ as well, $\log(1-p)<0$, and the ratio of two negative quantities is positive, so the lower bound on $B$ is a positive number. The ceiling rounds the real-valued lower bound upward. The limiting cases confirm the bound: if $q=0$, no finite positive budget can hit; if $q=1$, one draw already suffices since $1-(1-q)^B=1\ge p$ for any $B\ge 1$. This formula explains why random search can be effective when the good region is not too small. It also shows the limits: if $q$ is tiny, the required budget becomes enormous.

The coverage formula can be read as a budget-planning curve. The next cells evaluate $1-(1-q)^B$ over a range of budgets for three per-draw hit probabilities, then plot the resulting hit probabilities.