We fit a model on our data and get \(\hat{\beta}\)

Question: How well will this model predict a new observation \(y_{\text{new}}\) at \(x_{\text{new}}\)?

Fitted values: Predictions for observations in our dataset

\[\hat{y}_i = x_i^T\hat{\beta}\]

Out-of-sample predictions: Predictions for new observations

\[\hat{y}_{\text{new}} = x_{\text{new}}^T\hat{\beta}\]

Goal: Minimize prediction error on new data, not just training data!

Training error (in-sample): \[\text{MSE}_{\text{train}} = \frac{1}{n}\sum_{i=1}^n (y_i - \hat{y}_i)^2\]

Problem: Always decreases when we add predictors, even useless ones!

Solution: Estimate error on data the model hasn’t seen

Idea: For each observation \(i\):

1. Remove observation \(i\) from the dataset
2. Fit model on remaining \(n-1\) observations
   
3. Predict the removed observation
4. Calculate prediction error

\[\text{CV} = \frac{1}{n}\sum_{i=1}^n (y_i - \hat{y}_{-i})^2\]

where \(\hat{y}_{-i}\) is the prediction for \(y_i\) using all data except observation \(i\)

Problem: LOOCV requires fitting the model \(n\) times

For \(n = 1000\) observations, we need to fit 1000 models!

Good news: There’s a shortcut using the hat matrix!

Remember that fitted values are: \[\hat{y} = X\hat{\beta} = X(X^TX)^{-1}X^Ty = Hy\]

\[H = X(X^TX)^{-1}X^T\]

is called the hat matrix because it puts a “hat” on \(y\)

The diagonal elements \(h_{ii}\) are called leverages

\[\hat{y}_i = \sum_{j=1}^n h_{ij}y_j = h_{ii}y_i + \sum_{j \neq i}h_{ij}y_j\]

Interpretation: \(h_{ii}\) measures how much \(y_i\) influences its own prediction

Fun fact: We don’t need to refit the model \(n\) times!

\[\text{CV} = \frac{1}{n}\sum_{i=1}^n \left(\frac{y_i - \hat{y}_i}{1 - h_{ii}}\right)^2\]

This is exact LOOCV computed from a single model fit!

Intuition: The term \(\frac{1}{1-h_{ii}}\) inflates residuals for high-leverage points

When we remove observation \(i\), the prediction error is:

\[y_i - \hat{y}_{-i} = \frac{y_i - \hat{y}_i}{1 - h_{ii}}\]

Why? The fitted value \(\hat{y}_i\) uses \(y_i\) with weight \(h_{ii}\)

When we remove \(y_i\), we need to account for its influence

Idea: Penalize model complexity

\[\text{AIC} = n\log(\text{MSE}) + 2p+\underbrace{n+n\log(2\pi)}_{C}\]

where \(p\) is the number of parameters (including intercept)

Interpretation:

• First term: fit to the data (want small MSE)
  
• Second term: penalty for complexity (want small \(p\))

Lower AIC is better

Fun fact: AIC is asymptotically equivalent to LOOCV!

Implication: Choosing the model with smallest AIC is approximately the same as choosing the model with smallest LOOCV error

Notice: AIC and LOOCV select the same model (degree 2)!

Problem with \(R^2\): Always increases when we add predictors

\[R^2 = 1 - \frac{\text{SSE}}{\text{TSS}}\]

Solution: Penalize for number of predictors

\[R^2_{\text{adj}} = 1 - \frac{\text{SSE}/(n-p)}{\text{TSS}/(n-1)}\]

Higher is better (unlike AIC where lower is better)

Adjusted R² uses unbiased estimates of variance:

• \(\text{SSE}/(n-p)\) is unbiased for \(\sigma^2\)
  
• \(\text{TSS}/(n-1)\) is unbiased for variance of \(y\)
  

Adding a useless predictor:

• Decreases SSE slightly
  
• Increases \(p\) by 1
  
• Denominator \(n-p\) gets smaller
  
• Net effect can decrease \(R^2_{\text{adj}}\)

Another criterion for model selection:

\[C_p = \frac{\text{SSE}}{\hat{\sigma}^2} - n + 2p\]

where \(\hat{\sigma}^2\) is estimated from the largest model considered

Interpretation: Estimates prediction error, similar to AIC

So far: point prediction \(\hat{y}_{\text{new}} = x_{\text{new}}^T\hat{\beta}\)

Problem: This gives no sense of uncertainty!

Solution: Construct an interval that contains \(y_{\text{new}}\) with probability \(1-\alpha\)

Confidence interval for \(E[y_{\text{new}}]\):

• Uncertainty in estimating the mean response
  
• What is the average \(y\) at this \(x\)?

Prediction interval for \(y_{\text{new}}\):

• Uncertainty in predicting a single new observation
  
• What will the next \(y\) be at this \(x\)?

Prediction intervals are always wider!

For confidence intervals: \[\text{Var}(\hat{y}_{\text{new}}) = \sigma^2 x_{\text{new}}^T(X^TX)^{-1}x_{\text{new}}\]

For prediction intervals: \[\text{Var}(y_{\text{new}} - \hat{y}_{\text{new}}) = \sigma^2\left(1 + x_{\text{new}}^T(X^TX)^{-1}x_{\text{new}}\right)\]

The “+1” accounts for the inherent randomness of \(y_{\text{new}}\)

\[\hat{y}_{\text{new}} \pm t_{\alpha/2, n-p} \cdot \hat{\sigma}\sqrt{1 + x_{\text{new}}^T(X^TX)^{-1}x_{\text{new}}}\]

where:

• \(t_{\alpha/2, n-p}\) is the t-distribution critical value
  
• \(\hat{\sigma}\) is the estimated standard deviation of residuals

What happens when predictors are highly correlated?

• OLS estimates \(\hat{\beta}\) become unstable
• Small changes in data → large changes in coefficients
• Variance of \(\hat{\beta}\) is large
• Predictions can be unreliable

Quantifies multicollinearity:

\[\text{VIF}_j = \frac{1}{1 - R^2_j}\]

where \(R^2_j\) is from regressing \(x_j\) on all other predictors

Trade bias for reduced variance

Instead of minimizing SSE: \(\min_{\beta} \sum_{i=1}^n (y_i - x_i^T\beta)^2\)

Ridge regression adds a penalty: \(\min_{\beta} \sum_{i=1}^n (y_i - x_i^T\beta)^2 + \lambda\sum_{j=1}^p \beta_j^2\)

\(\lambda \geq 0\) is the tuning parameter:

• \(\lambda = 0\): ordinary least squares
  
• \(\lambda > 0\): shrinks coefficients toward zero
  
• \(\lambda \to \infty\): all coefficients → 0

Large coefficients indicate:

• Model is fitting noise
  
• Overly sensitive to small changes
  
• High variance in predictions

Penalizing \(\sum \beta_j^2\) encourages:

• Smaller, more stable coefficients
  
• Reduced variance
  
• Better prediction (despite bias)

The ridge estimator has a closed form: \(\hat{\beta}_{\text{ridge}} = (X^TX + \lambda I)^{-1}X^Ty\)

Compare to OLS: \(\hat{\beta}_{\text{OLS}} = (X^TX)^{-1}X^Ty\)

Key difference: Adding \(\lambda I\) to \(X^TX\) ensures the matrix is invertible even when \(X^TX\) is singular!

Problem: Penalty depends on scale of predictors

If \(x_1\) is in dollars and \(x_2\) is in thousands of dollars, their coefficients are penalized differently

Solution: Standardize predictors before applying ridge regression \(\tilde{x}_{ij} = \frac{x_{ij} - \bar{x}_j}{s_j}\)

Note: Usually don’t penalize the intercept

How do we choose \(\lambda\)?

Use cross-validation to find \(\lambda\) that minimizes prediction error!

LOOCV is a standard approach for selecting \(\lambda\)

Good news: The hat matrix trick works for ridge too!

For ridge regression with penalty \(\lambda\), the hat matrix is:

\(H_{\lambda} = X(X^TX + \lambda I)^{-1}X^T\)

LOOCV formula:

\(\text{CV}(\lambda) = \frac{1}{n}\sum_{i=1}^n \left(\frac{y_i - \hat{y}_i(\lambda)}{1 - h_{ii}(\lambda)}\right)^2\)

where \(h_{ii}(\lambda)\) is the \(i\)-th diagonal of \(H_{\lambda}\)

Just like OLS, but now everything depends on \(\lambda\)!

LOOCV advantages: - Exact formula using hat matrix - No randomness in the split - Maximum use of training data

K-fold CV (what cv.glmnet uses): - Computationally faster for large \(n\) - Can be more stable (less variance) - Default is 10-fold

In practice: Both work well, k-fold is more common

Ridge regression is most useful when:

1. Predictors are highly correlated
2. Number of predictors is large relative to sample size
3. You care more about prediction than interpretation
4. OLS estimates are unstable

Question: How many parameters does ridge regression use?

Naïve answer: Still \(p\) predictors… but they’re shrunk!

Correct answer: Effective degrees of freedom

\(\text{df}(\lambda) = \text{tr}(H_{\lambda})\)

where \(H_{\lambda} = X(X^TX + \lambda I)^{-1}X^T\) is the ridge hat matrix

Properties:

• \(\lambda = 0\): \(\text{df} = p\) (OLS)
  
• \(\lambda \to \infty\): \(\text{df} \to 0\) (no predictors).
• Larger \(\lambda\) → fewer effective parameters

Ridge solution: Point where SSE contour touches the constraint circle