{ "cells": [ { "cell_type": "raw", "metadata": {}, "source": [ "---\n", "title: 'Model selection and regularization'\n", "author: \"Sergio Bacallado, Jonathan Taylor\"\n", "subtitle: \"[web.stanford.edu/class/stats202](http://web.stanford.edu/class/stats202)\"\n", "date: \"Autumn 2020\"\n", "output:\n", " slidy_presentation:\n", " css: styles.css\n", "---" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's recap what we have learned about regression methods\n", "so far:\n", "\n", "- In linear regression, adding predictors always decreases the training error or RSS.\n", "\n", "- However, adding predictors does not necessarily improve the test error.\n", "\n", "- Selecting significant predictors is hard when $n$ is not much larger than $p$. \n", "\n", "- When $n\n", "\n", "\n", "\n", "## Best subset with `regsubsets`" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "library(ISLR) # where Credit is stored\n", "library(leaps) # where regsubsets is found\n", "summary(regsubsets(Balance ~ ., data=Credit))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "- Best model with 4 variables includes: `Cards, Income, Student, Limit`.\n", "\n", "## Choosing $k$ \n", "\n", "Naturally, $\\text{RSS}$ and \n", "$$\n", "R^2 = 1-\\frac{\\text{RSS}}{\\text{TSS}}\n", "$$\n", "improve as we increase $k$.\n", "\n", "To optimize $k$, we want to minimize the test error, not the training error. \n", "\n", "We could use **cross-validation**, or alternative estimates of test error:\n", "\n", "1. **Akaike Information Criterion (AIC)** (closely related to Mallow's $C_p$) given an estimate of the irreducible error $\\hat{\\sigma^2}$ :\n", "\n", "$$\\frac{1}{n}(\\text{RSS}+2k\\hat\\sigma^2)$$\n", "\n", "2. **Bayesian Information Criterion (BIC):**\n", "\n", "$$\\frac{1}{n}(\\text{RSS}+\\log(n)\\hat\\sigma^2)$$\n", "\n", "3. **Adjusted $R^2$**:\n", "\n", "$$R^2_a = 1-\\frac{\\text{RSS}/(n-k-1)}{\\text{TSS}/(n-1)}$$\n", "\n", "How do they compare to cross validation:\n", "\n", "- They are much less expensive to compute.\n", "\n", "- They are motivated by asymptotic arguments and rely on model assumptions (eg. normality of the errors).\n", "\n", "- Equivalent concepts for other models (e.g. logistic regression).\n", "\n", "## Example: best subset selection for the Credit dataset\n", "\n", "
\n", "\n", "
\n", "\n", "## Best subset selection for the Credit dataset\n", "\n", "
\n", "\n", "
\n", "\n", "*Recall:* In $K$-fold cross validation, we can estimate a standard error or accuracy for our test error estimate. Then, we can\n", "apply the *1SE rule*.\n", "\n", "# Stepwise selection methods\n", "\n", "Best subset selection has 2 problems:\n", "\n", "1. It is often very expensive computationally. We have to fit $2^p$ models!\n", "\n", "2. If for a fixed $k$, there are too many possibilities, we increase our chances of overfitting. The model selected has *high variance*." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In order to mitigate these problems, we can restrict our search space for the best model.\n", "This reduces the variance of the selected model at the expense of an increase in bias.\n", "\n", "## Forward selection\n", "\n", "
    \n", "
  1. Let ${\\cal M}_0$ denote the *null* model, which contains no predictors.\n", "
  2. For $k=0, \\dots, p-1$:\n", "
      \n", "
    1. Consider all $p-k$ models that augment the predictors\n", "in ${\\cal M}_k$ with one additional predictor.\n", "
    2. Choose the best among these $p-k$ models, and call it\n", "${\\cal M}_{k+1}$. Here *best* is defined as having smallest RSS or $R^2$.\n", "
    \n", "
  3. Select a single best model from among ${\\cal M}_0, \\dots, {\\cal M}_p$ using\n", "cross-validated prediction error, AIC, BIC or adjusted $R^2$.\n", "
\n", "\n", "## Forward selection in using `step`" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "M.full = lm(Balance ~ ., data=Credit)\n", "M.null = lm(Balance ~ 1, data=Credit)\n", "step(M.null, scope=list(upper=M.full), direction='forward', trace=0)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "- First four variables are `Rating, Income, Student, Limit`.\n", "\n", "- **Different from best subsets of size 4.**\n", "\n", "## Backward selection\n", "\n", "
    \n", "
  1. Let ${\\cal M}_p$ denote the *full* model, which contains all $p$ predictors\n", "
  2. For $k=p,p-1,\\dots,1$:\n", "
      \n", "
    1. Consider all $k$ models that contain\n", "all but one of the predictors in ${\\cal M}_k$ for a total of $k-1$ predictors.\n", "
    2. Choose the best among these $k$ models, and call it\n", "${\\cal M}_{k-1}$. Here *best* is defined as having smallest RSS or $R^2$.\n", "
    \n", "
  3. Select a single best model from among ${\\cal M}_0, \\dots, {\\cal M}_p$ using\n", "cross-validated prediction error, AIC, BIC or adjusted $R^2$.\n", "
" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Backward selection" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "step(M.full, scope=list(lower=M.null), direction='backward', trace=0)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Forward selection with BIC" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "step(M.null, scope=list(upper=M.full), direction='forward', trace=0, k=log(nrow(Credit))) # k defaults to 2, i.e AIC" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Forward vs. backward selection\n", "\n", "\n", "\n", "## Forward vs. backward selection" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "n = 5000\n", "X1 = rnorm(n)\n", "X2 = rnorm(n)\n", "X3 = X1 + 3 * X2 \n", "Y = X1 + 2 * X2 + rnorm(n)\n", "D = data.frame(X1, X2, X3, Y)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Forward" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "step(lm(Y ~ 1, data=D), list(upper=~ X1 + X2 + X3), direction='forward')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Backward" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "step(lm(Y ~ X1 + X2 + X3, data=D), list(lower=~ 1), direction='backward')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Other stepwise selection strategies\n", "\n", "- **Mixed stepwise selection (`direction='both'`):** Do forward selection, but at every step, remove any variables that are no longer *necessary*.\n", "\n", "- **Forward stagewise selection:** Roughly speaking, don't add in the variable *fully* at each step...\n", "\n", "- $\\dots$\n", "\n", "# Shrinkage methods\n", "\n", "**A mainstay of modern statistics!**\n", "\n", "The idea is to perform a linear regression, while *regularizing* or *shrinking* the coefficients $\\hat \\beta$ toward 0. \n", "\n", "## Why would shrunk coefficients be better?\n", "\n", "- This introduces *bias*, but may significantly decrease the *variance* of the estimates. If the latter effect is larger, this would decrease the test error.\n", "\n", "- Extreme example: set $\\hat{\\beta}$ to 0 -- variance is 0!\n", "\n", "- There are Bayesian motivations to do this: priors concentrated near 0 tend to shrink the parameters'\n", "posterior distribution. \n", "\n", "## Ridge regression\n", "\n", "Ridge regression solves the following optimization:\n", "\n", "$$\\min_{\\beta} \\;\\;\\; \\color{Blue}{ \\sum_{i=1}^n \\left(y_i -\\beta_0 -\\sum_{j=1}^p \\beta_jx_{i,j}\\right)^2} + \\color{Red}{\\lambda \\sum_{j=1}^p \\beta_j^2}$$\n", "\n", "- The RSS of the model at $\\beta$.\n", "\n", "- The squared $\\ell_2$ norm of $\\beta$, or $\\|\\beta\\|_2^2.$\n", "\n", "- The parameter $\\lambda$ is a tuning parameter. It modulates the importance of fit vs. shrinkage.\n", "\n", "- We find an estimate $\\hat\\beta^R_\\lambda$ for many values of $\\lambda$ and then choose it by cross-validation.\n", "Fortunately, this is no more expensive than running a least-squares regression.\n", "\n", "## Ridge regression\n", "\n", "\n", "\n", "## Ridge regression of `default` in the `Default` dataset\n", "\n", "
\n", "\n", "
\n", "\n", "## Ridge regression\n", "\n", "
\n", "\n", "
\n", "\n", "- In a simulation study, we compute squared bias, variance, and test error as a function of $\\lambda$.\n", "\n", "- In practice, cross validation would yield an estimate of the test error but bias and variance are unobservable.\n", "\n", "## Selecting $\\lambda$ by cross-validation for ridge regression\n", "\n", "
\n", "\n", "
\n", "\n", "## Lasso regression\n", "\n", "Lasso regression solves the following optimization:\n", "\n", "$$\\min_{\\beta} \\;\\;\\; \\color{Blue}{ \\sum_{i=1}^n \\left(y_i -\\beta_0 -\\sum_{j=1}^p \\beta_jx_{i,j}\\right)^2} + \\color{Red}{\\lambda \\sum_{j=1}^p |\\beta_j|}$$\n", "\n", "- RSS of the model at $\\beta$.\n", "\n", "- $\\ell_1$ norm of $\\beta$, or $\\|\\beta\\|_1.$\n", "\n", "- The parameter $\\lambda$ is a tuning parameter. It modulates the importance of fit vs. shrinkage.\n", "\n", "### Why would we use the Lasso instead of Ridge regression?\n", "\n", "- Ridge regression shrinks all the coefficients to a non-zero value.\n", "\n", "- The Lasso shrinks some of the coefficients all the way to zero.\n", "\n", "- A **convex** alternative (relaxation) to best subset selection (and its approximation stepwise selection)!\n", "\n", "## Ridge regression of `default` in the `Default` dataset\n", "\n", "
\n", "\n", "
\n", "\n", "A lot of pesky small coefficients throughout the regularization path\n", "\n", "## Lasso regression of `default` in the `Default` dataset\n", "\n", "
\n", "\n", "
\n", "\n", "Those coefficients are shrunk to zero\n", "\n", "## An alternative formulation for regularization\n", "\n", "