---
title: Optimism and Over-Fitting
author: 36-465/665, Spring 2021
date: 11 February 2021
output: slidy_presentation
---

\[
\newcommand{\Prob}[1]{\mathbb{P}\left( #1 \right)}
\newcommand{\Expect}[1]{\mathbb{E}\left[ #1 \right]}
\newcommand{\Var}[1]{\mathrm{Var}\left[ #1 \right]}
\newcommand{\Cov}[1]{\mathrm{Cov}\left[ #1 \right]}
\newcommand{\Risk}{r}
\newcommand{\EmpRisk}{\hat{r}}
\newcommand{\Loss}{\ell}
\newcommand{\OptimalStrategy}{\sigma}
\DeclareMathOperator*{\argmin}{argmin}
\newcommand{\ModelClass}{S}
\newcommand{\OptimalModel}{s^*}
\DeclareMathOperator{\tr}{tr}
\]

```{r, include=FALSE}
# General set-up options
library(knitr)
opts_chunk$set(size="small", background="white", highlight=FALSE,
               cache=TRUE, autodep=TRUE,
               tidy=TRUE, warning=FALSE, message=FALSE,
               echo=FALSE)
set.seed(2020-03-20)
```


## Housekeeping

- HW 1 is due this evening at 6 pm (Pittsburgh time)
- Also be sure to do HW 0 about course policies (on Canvas, auto-graded)

## In our last episode

> - True risk $r(s) \equiv \Expect{\Loss(Y, s(X))}$
> - Empirical risk $\EmpRisk(s) \equiv n^{-1}\sum_{i=1}^{n}{\Loss(y_i, s(x_i))}$
> - Empirical risk minimization: find $\hat{s}$ that minimizes $\hat{r}$ and use it to predict
> - With IID data, $\EmpRisk(s)$ is an unbiased estimate of $\Risk(s)$ for any _fixed_ $s$
> - and $\EmpRisk(s) \rightarrow \Risk(s)$ as $n\rightarrow\infty$ (law of large numbers, again for fixed $s$)
> - and $\Var{\EmpRisk(s)} = O(1/n)$ (again for fixed $s$)
> - "Parameterize" strategies with vector $\theta$; then $\hat{\theta} \rightarrow \theta^*$, and usually $\hat{\theta} = \theta^* + O(1/\sqrt{n})$
> - But parameters don't matter, predictions do...

## Reminders about "the usual asymptotics"

- Upshots of the [Taylor expansion song and dance from last time](http://www.stat.cmu.edu/~cshalizi/sml/21/lectures/03/lecture-03.html#(14)):

\begin{eqnarray}
\hat{\theta} & \rightarrow & \theta^* \\
\hat{\theta} & \approx & \theta^* - \mathbf{k}^{-1} \nabla \EmpRisk(\theta^*)\\
\mathbf{k} & \equiv & \nabla \nabla \Risk(\theta^*)\\
\Var{\hat{\theta}} & \approx & n^{-1} \mathbf{k}^{-1} \mathbf{j} \mathbf{k}^{-1}\\
\mathbf{j} & \equiv & \Var{\nabla \Loss(Y, s(X;\theta^*))}\\
\hat{\theta} & \rightsquigarrow & \mathcal{N}(\theta^*, n^{-1} \mathbf{k}^{-1} \mathbf{j} \mathbf{k}^{-1})
\end{eqnarray}

## How do we find those magic matrices?

- Approximations/estimates, since $\EmpRisk \rightarrow \Risk$ and $\hat{\theta} \rightarrow \theta^*$:
\begin{eqnarray}
\mathbf{k} & \equiv & \nabla\nabla \Risk(\theta^*)\\
& \approx & \nabla \nabla \Risk(\hat{\theta})\\
& \approx & \nabla \nabla \EmpRisk(\hat{\theta})\\
\mathbf{j} & \equiv & \Var{\nabla \Loss(Y, s(X;\theta^*))}\\
& \approx & \Var{\nabla \Loss(Y, s(X;\hat{\theta}))}\\
& \approx & \frac{1}{n}\sum_{i=1}^{n}{\left( \nabla \Loss(y_i, s(x_i;\hat{\theta}))\right) \left( \nabla \Loss(y_i, s(x_i;\hat{\theta}))\right)^T}
\end{eqnarray}
    + (last equation expression $=$ sample variance-covariance matrix for $\nabla \Loss$)
- Just need to compute $\nabla \Loss$ and $\nabla \nabla \Loss$
    + Often we did that to find $\hat{\theta}$ in the first place
- You'll work through an example in HW 2


## But what about the predictions?

- $\theta^*$ is the best parameter value, so best predictions are $s(x;\theta^*) = s^*(x)$
    + "Best" = risk-minimizing = expected loss minimizing
- $\hat{\theta}$ estimates $\theta^*$, so $s(x;\hat{\theta}) = \hat{s}(x)$ estimates the best prediction
- Uncertainty about $\theta^*$ $\Rightarrow$ uncertainty about best predictions
- Variance of $\hat{\theta}$ $\Rightarrow$ variance for best predictions
    + "Propagation of error", "uncertainty propagation", "the delta method"
	+ Another problem in HW 2
	+ More fun with $\nabla$s and matrices

## But what about the predictions? (2)

- To preview the homework a little, use a Taylor expansion:
\begin{eqnarray}
s(x;\hat{\theta}) & \approx & s(x;\theta^*) + \nabla s(x;\theta^*) (\hat{\theta} - \theta^*)\\
& = & s(x;\theta^*) +  \nabla s(x;\theta^*) O(1/\sqrt{n})\\
& = & s(x;\theta^*) + O(1/\sqrt{n}) \rightarrow s(x;\theta^*)\\
\end{eqnarray}

> - We'll be more precise in the homework, rather than just $O(1/\sqrt{n})$


## What about the risk?

- Taylor expand again:
\begin{eqnarray}
\Risk(\hat{\theta}) & \approx & \Risk(\theta^*) + \frac{1}{2} (\hat{\theta} - \theta^*) \cdot \mathbf{k} (\hat{\theta} - \theta^*)\\
& = & \Risk(\theta^*) + O(1/\sqrt{n})O(1/\sqrt{n})\\
& = & \Risk(\theta^*) + O(1/n)
\end{eqnarray}
- Remember that in [Lecture 2](http://www.stat.cmu.edu/~cshalizi/sml/21/lectures/02/lecture-02.html#(11)) we said
\[
\Risk(s) = \Risk_0 + (\Risk(s^*) - \Risk_0) + (\Risk(s) - \Risk(s^*)) = \text{minimum risk} + \text{approximation error} + \text{estimation error}
\]
- Generically, estimation error for ERM is $O(1/n)$
    + Details later today


## Being more precise about the risk

- For any _fixed_ $s$,
\[
\Expect{\EmpRisk(s)} = \Risk(s)
\]
- Define $\gamma(s) \equiv \EmpRisk(s) - \Risk(s)$, the **deviation** or **fluctuation** at $s$
- $\Expect{\gamma(s)} = 0$ for any fixed $s$, and $\Var{\gamma(s)} = O(1/n)$
- But what about $\gamma(\hat{s})$?  $\hat{s}$ is _not_ a fixed strategy!!!

## Empirical risk minimization is optimistic

> - Intuition: $\hat{s}$ picked to fit _this_ data, so it partly adapts to real patterns that will show up in new data, and partly to the accidents of the training data
> - Math:
\[
\hat{s} = \argmin_{s \in S}{\left( \Risk(s) + \gamma(s) \right)}
\]
> - Picking $s$ to minimize $\EmpRisk(s)$ is _partly_ about finding the strategy with the smallest true risk, minimizing $\Risk(s)$
> - Picking $s$ to minimize $\EmpRisk(s)$ is _partly_ about finding the _luckiest_ strategy, one with very negative $\gamma(s)$
>     + Remember $\Expect{\gamma(s)} = 0$ for each $s$
> - Implication: $\Expect{\gamma(\hat{s})} < 0$
> - Implication:
\begin{eqnarray}
\Expect{\EmpRisk(\hat{s})} & = & \Expect{\Risk(\hat{s})} + \Expect{\gamma(\hat{s})}\\
& < & \Expect{\Risk(\hat{s})}
\end{eqnarray}
> - **The strategy we pick by minimizing the in-sample loss will do worse on new data (on average)**


## A toy example (1)

```{r}
# Demo data
y <- rnorm(n=20, mean=7, sd=1)
```

- $Y \sim \mathcal{N}(7, 1)$, $n=20$ data points (tick marks on the axis):

```{r demo.of.erm}
plot(0, type="n", xlim=range(y), ylim=c(0, 5), xlab=expression(theta),
     ylab=expression(r(theta)))
rug(x=y, side=1)
```

## A toy example (2)

- Use squared error, with $\ModelClass =$ all constants, so we just want to predict the expected value of $Y$
- True risk when predicting $\theta$ is $1+(\theta-7)^2$ (why?)
  + What if $Y$ is Gaussian but with a different variance?  What if it was non-Gaussian but with the same expectation and variance?
- True risk is minimized at the expected value:

```{r}
<<demo.of.erm>>
curve(1+(x-7)^2, add=TRUE, col="grey")
legend("bottomright", legend=expression(r(theta)), col="grey", lty=1)
```


## A toy example (3)

- _Empirical_ risk is minimized at the sample mean:

```{r}
<<demo.of.erm>>
curve(1+(x-7)^2, add=TRUE, col="grey")
# curve() is a bit fussy about needing "vectorized" functions
empirical.risk <- function(m) { mean((y-m)^2) }
curve(Vectorize(empirical.risk)(x), add=TRUE, col="black")
legend("bottomright", legend=c(expression(r(theta)),
                               expression(hat(r)(theta))),
       col=c("grey", "black"), lty=1)
```


## A toy example (4)

- The difference between the two curves is the risk deviation $\gamma(\theta)$:
```{r}
curve(Vectorize(empirical.risk)(x) - (1+(x-7)^2),
      from=min(y), to=max(y),
      xlab=expression(theta), ylab=expression(gamma(theta)))
abline(v=mean(y), lty="dotted")
abline(v=7, lty="dotted", col="grey")
```

- $\gamma(\theta)$ is really a random _function_ (a **stochastic process**), and this is one draw from its distribution (one **realization** of the process)


## A toy example (5)

- Repeat the simulation many times, to evaluate $\gamma(\theta)$ for any fixed $\theta$
    + At each $\theta$, get a Gaussian distribution _centered at 0_ (why?)


```{r}
# Risk deviations at an arbitrary _fixed_ parameter value
gamma.fixed.m <- function(m, b=2000) {
    replicate(b, mean((rnorm(n=20, mean=7, sd=1) - m)^2) - (1+(7-m)^2))
}    
```
```{r}
# Risk deviations at the optimal parameter value
gamma.7 <- gamma.fixed.m(7)
hist(gamma.7, xlab=expression(gamma(7)),
     main="Risk deviations at the optimal parameter value (theta=7)",
     sub=paste("Average deviation over 2000 simulations of 20 draws each=", signif(mean(gamma.7), 2)))
```
```{r}
# Risk deviations at an arbitrary, non-optimal, fixed parameter vale
gamma.4.2 <- gamma.fixed.m(4.2)
hist(gamma.4.2, xlab=expression(gamma(4.2)),
     main="Risk deviations at a fixed, non-optimal parameter value (theta=4.2)",
     sub=paste("Average deviation over 2000 simulations of 20 draws each=", signif(mean(gamma.4.2), 2)))
```


## A toy example (6)

- Things look different for the $\hat{\theta}$ chosen by empirical risk minimization:

```{r}
# Risk deviations for the parameter that minimizes in-sample risk
gamma.erm <- function(b=2000) {
    gamma.erm.once <- function() {
        y.new <- rnorm(n=20, mean=7, sd=1)
        hat.m <- mean(y.new)
        return(mean((y.new-hat.m)^2) - (1+(7-hat.m)^2))
    }
    replicate(b, gamma.erm.once())
}

gamma.adaptive <- gamma.erm()
hist(gamma.adaptive, xlab=expression(gamma(hat(theta))),
     main="Risk deviations for the parameter value that minimizes empirical risk",
     sub=paste("Average deviation over 2000 simulations of 20 draws each=", signif(mean(gamma.adaptive),2)))
```


## Estimating the true risk from the empirical risk

- We want to know $\Risk(\hat{s})$, and we know $\EmpRisk(\hat{s})$
    + Equivalently, $\Risk(\hat{\theta})$ and $\EmpRisk(\hat{\theta})$
- We saw a little while ago that
\[
\Risk(\hat{\theta})  \approx \Risk(\theta^*) + \frac{1}{2} (\hat{\theta} - \theta^*) \cdot \mathbf{k} (\hat{\theta} - \theta^*)\\
\]
- We don't know $\Risk(\theta^*)$ but $\theta^*$ _is_ fixed so $\Risk(\theta^*) = \Expect{\EmpRisk(\theta^*)}$
- We don't know $\EmpRisk(\theta^*)$ but we can Taylor expand:
\begin{eqnarray}
\EmpRisk(\theta^*) & \approx & \EmpRisk(\hat{\theta}) + \frac{1}{2} (\theta^* - \hat{\theta}) \cdot \left( \nabla \nabla \EmpRisk(\hat{\theta}) (\theta^* - \hat{\theta})\right)\\
& = & \EmpRisk(\hat{\theta}) + \frac{1}{2} (\hat{\theta} - \theta^*) \cdot \left( \nabla \nabla \EmpRisk(\hat{\theta}) (\hat{\theta} - \theta^*)\right)\\
& \approx & \EmpRisk(\hat{\theta}) + \frac{1}{2} (\hat{\theta} - \theta^*) \cdot \left( \nabla \nabla \Risk(\hat{\theta}) (\hat{\theta} - \theta^*)\right)\\
& = & \EmpRisk(\hat{\theta}) + \frac{1}{2} (\hat{\theta} - \theta^*) \cdot \left( \mathbf{k} (\hat{\theta} - \theta^*)\right)\\
\end{eqnarray}

## Estimating the true risk from the empirical risk (2)

- Put the two Taylor series together
\[
\Risk(\hat{\theta}) \approx \EmpRisk(\hat{\theta}) + (\hat{\theta} - \theta^*) \cdot \mathbf{k} (\hat{\theta} - \theta^*)
\]
- Now (see back-up)
\[
\Expect{\Risk(\hat{\theta})} \approx \Expect{\EmpRisk(\hat{\theta})} + n^{-1}\tr{\left(\mathbf{j}\mathbf{k}^{-1}\right)}
\]
    + Remember $\tr{\mathbf{a}} \equiv \sum_{i}{a_{ii}}=$ sum of the diagonal entries of $\mathbf{a}$ $=$ sum of the eigenvalues of $\mathbf{a}$
- So we can _estimate_
\[
\Risk(\hat{\theta}) \approx \EmpRisk(\hat{\theta}) + n^{-1}\tr{\mathbf{j}\mathbf{k}^{-1}}
\]


## Estimating the true risk from the empirical risk (3)

\[
\Risk(\hat{\theta}) \approx \EmpRisk(\hat{\theta}) + n^{-1}\tr{\mathbf{j}\mathbf{k}^{-1}}
\]

- An (asymptotically) unbiased estimate of the true risk based on the empirical risk
    + Still true even if we're more careful in our math
- The extra term is always $> 0$
- The extra term measures how "flat" the minimum is ($\mathbf{k}^{-1}$) and how much noise there is in the function we're optimizing ($\mathbf{j}$)
- Lots of other formulas are special cases of this:
    + Akaike information criterion (AIC) for model selection with maximum likelihood
	    * If we use the log loss _and_ our model is right, _then_ $\tr{(\mathbf{j}\mathbf{k}^{-1})} =$ number of dimensions in $\theta$
	+ Mallows $C_p$ for linear regression, $\EmpRisk(\hat{\theta}) + 2\sigma^2 p/n$ (when we estimate $p$ coefficients and the true noise around the regression line has variance $\sigma^2$)
	+ Generalized Mallows $C_p$ for linear smoothers, $\EmpRisk(\hat{s}) + 2\sigma^2 \mathrm{df}(\hat{s})/n$
	    * $\mathrm{df}(\hat{s}) \equiv \sum_{i=1}^{\sigma^2}{\Cov{Y_i, s(X_i, \hat{\theta})}}$
- The special cases let us see the $\tr{(\mathbf{j}\mathbf{k}^{-1})}$ term is _something like_ "how flexible is the model?"
	




## The moral of all this math: approximation-estimation trade-off {.smaller}

- Remember how we broke up the risk in [Lecture 2](http://www.stat.cmu.edu/~cshalizi/sml/21/lectures/02/lecture-02.html#(11)):
\begin{eqnarray}
\Risk(\hat{s}) & = & \Risk_0 + (\Risk(\OptimalModel) - \Risk_0) + (\Risk(\hat{s}) - \Risk(\OptimalModel))\\
& = & \text{(optimal risk)} + \text{(approximation error)} + \text{(estimation error)}
\end{eqnarray}
- What happens if we have _two_ sets of strategies, $\ModelClass_1$ and $\ModelClass_2$, and $\ModelClass_1 \subset \ModelClass_2$?
- Usually two optimal models, $\OptimalModel_1 \in \ModelClass_1$ and $\OptimalModel_2 \in \ModelClass_2$, and two estimates, $\hat{s}_1$ and $\hat{s}_2$
- Because $\ModelClass_1 \subset \ModelClass_2$:
    1. $\EmpRisk(\hat{s}_1) \geq \EmpRisk(\hat{s}_2)$: empirical risk can only get better by optimizing over more strategies
	2. $\Risk(\OptimalModel_1) \geq \Risk(\OptimalModel_2)$: true risk can only get better by optimizing over more strategies
	3. $\max_{s\in\ModelClass_1}{|\gamma(s)|} \leq \max_{s\in\ModelClass_2}{|\gamma(s)|}$: the maximum deviation can only get bigger by searching over more strategies
- Thing (1) means ERM always prefers the bigger, more flexible model class
- Thing (2) means that bigger, more flexible model classes will have smaller approximation error
- Thing (3) means that bigger, more flexible model classes will have bigger estimation error

## Over-fitting

- **Over-fitting** = fitting a model that's bigger (more flexible, more powerful) than the one which will predict best (given your data)
- Over-fitting happens _because_ optimism increases with model size
    + Equivalently, _because_ estimation error grows with model size
- Some amount of optimism is built in to empirical risk minimization
    + We _asked_ the strategy to optimize performance on the training data, and it did so
- ERM is _very_ prone to over-fitting

## Avoiding over-fitting

> - Good practical advice: constrain the strategies using what you know (or believe...) about the situation
> - Good practical advice: keep it simple
>    + Tension: sometimes we know things _are_ complicated (e.g. language modeling)
> - Impose penalties, so we minimize $\EmpRisk(s) + g(s)$ for some $g(s)$ which is something like $\Expect{|\gamma(s)|}$
>    + Don't take that too literally; we'll come back to penalties later
> - Get tighter control on over-fitting


## Why the optimism is not the end of the story

\[
\Risk(\hat{\theta}) \approx \EmpRisk(\hat{\theta}) + n^{-1}\tr{\mathbf{j}\mathbf{k}^{-1}}
\]

- An (asymptotically) unbiased estimate of the true risk based on the empirical risk
- What about the variance?  What about the distribution?
- What about _non_-asymptotic results?


## The goal: generalization error bounds

- Suppose we could say: for any $n$ and $\alpha > 0$,
\[
\Prob{\Risk(\hat{\theta}) - \EmpRisk(\hat{\theta}) \geq g(n, \alpha)} \leq \alpha
\]
for some function $g(n,\alpha)$
    + Presumably $g(n,\alpha) \rightarrow 0$ as $n\rightarrow \infty$ for constant $\alpha$
- This'd be a _bound_ on the **generalization error**, a **generalization error bound**
    + Or: a one-sided confidence interval for $\Risk(\hat{\theta})$, confidence level $=1-\alpha$
- Generalization error bounds would control over-fitting
- Over the next few weeks we'll see how to find them:
    1. Bounding the deviation of a _single_ function from its expectation value
	2. Bounding the _maximum_ deviation over _infinitely many_ functions
	    * Remember the toy example
	3. Measuring model flexibility or complexity
	
	
## Summing up

- The empirical risk of a _fixed_ strategy is an unbiased estimate of its true risk
- The strategy $\hat{s}$ picked by ERM isn't fixed but random
- The empirical risk of $\hat{s}$ is a _negatively-biased_ or optimistic estimate of the true risk of $\hat{s}$
- Larger, more flexible models have smaller approximation error, but more optimism and so more over-fitting
- The usual asymptotics give us rough, asymptotic formulas for the optimism / generalization error
- We'll see how to get non-asymptotic generalization error bounds







## Backup: Expectations of quadratic forms, and the optimism formula

- We want
\[
\Expect{(\hat{\theta} - \theta^*) \cdot \left( \mathbf{k} (\hat{\theta} - \theta^*)\right)}
\]
- A **quadratic form** is an expression like $Z \cdot \mathbf{b} Z$, where $Z$ is a random vector and $\mathbf{b}$ is a _non_-random matrix
    + $\Expect{Z} = \mu$ (a vector) and $\Var{Z} = \Sigma$ (a matrix)
    + The quadratic form is a random scalar
- Fun math fact:
\[
\Expect{Z \cdot \mathbf{b} Z} = \tr{(\mathbf{b} \mathbf{\Sigma})} + \mu \cdot \mathbf{b} \mu
\]
    + It will be character-building to prove this without looking it up (hint: trace is cyclic)
- Applied here, $Z = \hat{\theta}-\theta^*$
    + $\Expect{Z} = \Expect{\hat{\theta}} - \theta^* \rightarrow 0$ (as $n\rightarrow\infty$)
	+ $\Var{Z} = n^{-1} \mathbf{k}^{-1}\mathbf{j}\mathbf{k}^{-1}$ (as $n\rightarrow\infty$)
- So
\[
\Expect{(\hat{\theta} - \theta^*) \cdot \left( \mathbf{k} (\hat{\theta} - \theta^*)\right)} \rightarrow \tr{(\mathbf{k} n^{-1} \mathbf{k}^{-1} \mathbf{j} \mathbf{k}^{-1})} = n^{-1} \tr{(\mathbf{j}\mathbf{k}^{-1})}
\]

