cross-validation-bootstrap.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Cross validation and bootstrap lab (R)\n",
    "\n",
    "From An Introduction to Statistical Learning, chapter 5 lab."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 51,
   "metadata": {},
   "outputs": [],
   "source": [
    "library(ISLR)\n",
    "library(boot)\n",
    "\n",
    "set.seed(1)\n",
    "train = sample(392, 196)\n",
    "\n",
    "#attach(Auto)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1 Validation set\n",
    "\n",
    "### 1.1 Linear fit"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "23.2660086465003"
      ],
      "text/latex": [
       "23.2660086465003"
      ],
      "text/markdown": [
       "23.2660086465003"
      ],
      "text/plain": [
       "[1] 23.26601"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "lm.fit = lm(mpg~horsepower, data=Auto, subset=train)\n",
    "mean((mpg - predict(lm.fit, Auto))[-train]^2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.2 Quadratic fit"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "18.7164594933828"
      ],
      "text/latex": [
       "18.7164594933828"
      ],
      "text/markdown": [
       "18.7164594933828"
      ],
      "text/plain": [
       "[1] 18.71646"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "lm.fit2 = lm(mpg~poly(horsepower, 2), data=Auto, subset=train)\n",
    "mean((mpg - predict(lm.fit2, Auto))[-train]^2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.3 Cubic fit"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "18.7940067973945"
      ],
      "text/latex": [
       "18.7940067973945"
      ],
      "text/markdown": [
       "18.7940067973945"
      ],
      "text/plain": [
       "[1] 18.79401"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "lm.fit3 = lm(mpg~poly(horsepower, 3), data=Auto, subset=train)\n",
    "mean((mpg - predict(lm.fit3, Auto))[-train]^2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2 Leave-one-out cross validation (LOOCV)\n",
    "\n",
    "### 2.1 Linear fit"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<style>\n",
       ".dl-inline {width: auto; margin:0; padding: 0}\n",
       ".dl-inline>dt, .dl-inline>dd {float: none; width: auto; display: inline-block}\n",
       ".dl-inline>dt::after {content: \":\\0020\"; padding-right: .5ex}\n",
       ".dl-inline>dt:not(:first-of-type) {padding-left: .5ex}\n",
       "</style><dl class=dl-inline><dt>(Intercept)</dt><dd>39.9358610211705</dd><dt>horsepower</dt><dd>-0.157844733353654</dd></dl>\n"
      ],
      "text/latex": [
       "\\begin{description*}\n",
       "\\item[(Intercept)] 39.9358610211705\n",
       "\\item[horsepower] -0.157844733353654\n",
       "\\end{description*}\n"
      ],
      "text/markdown": [
       "(Intercept)\n",
       ":   39.9358610211705horsepower\n",
       ":   -0.157844733353654\n",
       "\n"
      ],
      "text/plain": [
       "(Intercept)  horsepower \n",
       " 39.9358610  -0.1578447 "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "glm.fit = glm(mpg~horsepower, data=Auto)\n",
    "coef(glm.fit)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 56,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<style>\n",
       ".list-inline {list-style: none; margin:0; padding: 0}\n",
       ".list-inline>li {display: inline-block}\n",
       ".list-inline>li:not(:last-child)::after {content: \"\\00b7\"; padding: 0 .5ex}\n",
       "</style>\n",
       "<ol class=list-inline><li>19.0332138547041</li><li>19.0324172666451</li></ol>\n"
      ],
      "text/latex": [
       "\\begin{enumerate*}\n",
       "\\item 19.0332138547041\n",
       "\\item 19.0324172666451\n",
       "\\end{enumerate*}\n"
      ],
      "text/markdown": [
       "1. 19.0332138547041\n",
       "2. 19.0324172666451\n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "[1] 19.03321 19.03242"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "cv.err = cv.glm(Auto, glm.fit)\n",
    "cv.err$delta"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.2 Degree 1-5 fits"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<style>\n",
       ".list-inline {list-style: none; margin:0; padding: 0}\n",
       ".list-inline>li {display: inline-block}\n",
       ".list-inline>li:not(:last-child)::after {content: \"\\00b7\"; padding: 0 .5ex}\n",
       "</style>\n",
       "<ol class=list-inline><li>24.2315135179292</li><li>19.2482131244897</li><li>19.3349840640291</li><li>19.4244303104302</li><li>19.0332138547041</li></ol>\n"
      ],
      "text/latex": [
       "\\begin{enumerate*}\n",
       "\\item 24.2315135179292\n",
       "\\item 19.2482131244897\n",
       "\\item 19.3349840640291\n",
       "\\item 19.4244303104302\n",
       "\\item 19.0332138547041\n",
       "\\end{enumerate*}\n"
      ],
      "text/markdown": [
       "1. 24.2315135179292\n",
       "2. 19.2482131244897\n",
       "3. 19.3349840640291\n",
       "4. 19.4244303104302\n",
       "5. 19.0332138547041\n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "[1] 24.23151 19.24821 19.33498 19.42443 19.03321"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "cv.error = rep(0, 5)\n",
    "for (i in 1:5) {\n",
    "    glm.fit = glm(mpg~poly(horsepower, i), data=Auto)\n",
    "    cv.error[i] = cv.glm(Auto, glm.fit)$delta[1]\n",
    "}\n",
    "cv.error"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3 k-fold cross validation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<style>\n",
       ".list-inline {list-style: none; margin:0; padding: 0}\n",
       ".list-inline>li {display: inline-block}\n",
       ".list-inline>li:not(:last-child)::after {content: \"\\00b7\"; padding: 0 .5ex}\n",
       "</style>\n",
       "<ol class=list-inline><li>24.2720671232254</li><li>19.2690928085129</li><li>19.3480535605547</li><li>19.2949648229745</li><li>19.0319790002896</li><li>18.8978121056401</li><li>19.1206066690695</li><li>19.1466631054789</li><li>18.8701307442148</li><li>20.9552042280394</li></ol>\n"
      ],
      "text/latex": [
       "\\begin{enumerate*}\n",
       "\\item 24.2720671232254\n",
       "\\item 19.2690928085129\n",
       "\\item 19.3480535605547\n",
       "\\item 19.2949648229745\n",
       "\\item 19.0319790002896\n",
       "\\item 18.8978121056401\n",
       "\\item 19.1206066690695\n",
       "\\item 19.1466631054789\n",
       "\\item 18.8701307442148\n",
       "\\item 20.9552042280394\n",
       "\\end{enumerate*}\n"
      ],
      "text/markdown": [
       "1. 24.2720671232254\n",
       "2. 19.2690928085129\n",
       "3. 19.3480535605547\n",
       "4. 19.2949648229745\n",
       "5. 19.0319790002896\n",
       "6. 18.8978121056401\n",
       "7. 19.1206066690695\n",
       "8. 19.1466631054789\n",
       "9. 18.8701307442148\n",
       "10. 20.9552042280394\n",
       "\n",
       "\n"
      ],
      "text/plain": [
       " [1] 24.27207 19.26909 19.34805 19.29496 19.03198 18.89781 19.12061 19.14666\n",
       " [9] 18.87013 20.95520"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "set.seed(17)\n",
    "cv.error.10 = rep(0,10)\n",
    "for (i in  1:10) {\n",
    "    glm.fit = glm(mpg~poly(horsepower, i), data=Auto)\n",
    "    cv.error.10[i] = cv.glm(Auto, glm.fit, K=10)$delta[1]\n",
    "}\n",
    "cv.error.10"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4 The bootstrap\n",
    "\n",
    "### 4.1 Estimating the accuracy of a statistic (alpha)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 60,
   "metadata": {},
   "outputs": [],
   "source": [
    "alpha.fn = function(data, index) {\n",
    "    X = data$X[index]\n",
    "    Y = data$Y[index]\n",
    "    return (var(Y) - cov(X, Y)) / (var(X) + var(Y) - 2 * cov(X, Y))\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 61,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "0.681879176718024"
      ],
      "text/latex": [
       "0.681879176718024"
      ],
      "text/markdown": [
       "0.681879176718024"
      ],
      "text/plain": [
       "[1] 0.6818792"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "alpha.fn(Portfolio, 1:100)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 63,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "0.774831292847958"
      ],
      "text/latex": [
       "0.774831292847958"
      ],
      "text/markdown": [
       "0.774831292847958"
      ],
      "text/plain": [
       "[1] 0.7748313"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "set.seed(1)\n",
    "alpha.fn(Portfolio, sample(100, 100, replace=T))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 64,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "\n",
       "ORDINARY NONPARAMETRIC BOOTSTRAP\n",
       "\n",
       "\n",
       "Call:\n",
       "boot(data = Portfolio, statistic = alpha.fn, R = 1000)\n",
       "\n",
       "\n",
       "Bootstrap Statistics :\n",
       "     original      bias    std. error\n",
       "t1* 0.6818792 -0.01141145   0.1235559"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "boot(Portfolio, alpha.fn, R=1000)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 4.2 Estimating the accuracy of a linear regression model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 66,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<style>\n",
       ".dl-inline {width: auto; margin:0; padding: 0}\n",
       ".dl-inline>dt, .dl-inline>dd {float: none; width: auto; display: inline-block}\n",
       ".dl-inline>dt::after {content: \":\\0020\"; padding-right: .5ex}\n",
       ".dl-inline>dt:not(:first-of-type) {padding-left: .5ex}\n",
       "</style><dl class=dl-inline><dt>(Intercept)</dt><dd>39.9358610211705</dd><dt>horsepower</dt><dd>-0.157844733353654</dd></dl>\n"
      ],
      "text/latex": [
       "\\begin{description*}\n",
       "\\item[(Intercept)] 39.9358610211705\n",
       "\\item[horsepower] -0.157844733353654\n",
       "\\end{description*}\n"
      ],
      "text/markdown": [
       "(Intercept)\n",
       ":   39.9358610211705horsepower\n",
       ":   -0.157844733353654\n",
       "\n"
      ],
      "text/plain": [
       "(Intercept)  horsepower \n",
       " 39.9358610  -0.1578447 "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "boot.fn = function(data, index) {\n",
    "    return (coef(lm(mpg~horsepower, data=data, subset=index)))\n",
    "}\n",
    "boot.fn(Auto, 1:392)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 67,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<style>\n",
       ".dl-inline {width: auto; margin:0; padding: 0}\n",
       ".dl-inline>dt, .dl-inline>dd {float: none; width: auto; display: inline-block}\n",
       ".dl-inline>dt::after {content: \":\\0020\"; padding-right: .5ex}\n",
       ".dl-inline>dt:not(:first-of-type) {padding-left: .5ex}\n",
       "</style><dl class=dl-inline><dt>(Intercept)</dt><dd>40.3404516830189</dd><dt>horsepower</dt><dd>-0.163486837689938</dd></dl>\n"
      ],
      "text/latex": [
       "\\begin{description*}\n",
       "\\item[(Intercept)] 40.3404516830189\n",
       "\\item[horsepower] -0.163486837689938\n",
       "\\end{description*}\n"
      ],
      "text/markdown": [
       "(Intercept)\n",
       ":   40.3404516830189horsepower\n",
       ":   -0.163486837689938\n",
       "\n"
      ],
      "text/plain": [
       "(Intercept)  horsepower \n",
       " 40.3404517  -0.1634868 "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "set.seed(1)\n",
    "boot.fn(Auto, sample(392, 392, replace=T))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 68,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<style>\n",
       ".dl-inline {width: auto; margin:0; padding: 0}\n",
       ".dl-inline>dt, .dl-inline>dd {float: none; width: auto; display: inline-block}\n",
       ".dl-inline>dt::after {content: \":\\0020\"; padding-right: .5ex}\n",
       ".dl-inline>dt:not(:first-of-type) {padding-left: .5ex}\n",
       "</style><dl class=dl-inline><dt>(Intercept)</dt><dd>40.1186906449022</dd><dt>horsepower</dt><dd>-0.157706320543503</dd></dl>\n"
      ],
      "text/latex": [
       "\\begin{description*}\n",
       "\\item[(Intercept)] 40.1186906449022\n",
       "\\item[horsepower] -0.157706320543503\n",
       "\\end{description*}\n"
      ],
      "text/markdown": [
       "(Intercept)\n",
       ":   40.1186906449022horsepower\n",
       ":   -0.157706320543503\n",
       "\n"
      ],
      "text/plain": [
       "(Intercept)  horsepower \n",
       " 40.1186906  -0.1577063 "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "boot.fn(Auto, sample(392, 392, replace=T))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 69,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "\n",
       "ORDINARY NONPARAMETRIC BOOTSTRAP\n",
       "\n",
       "\n",
       "Call:\n",
       "boot(data = Auto, statistic = boot.fn, R = 1000)\n",
       "\n",
       "\n",
       "Bootstrap Statistics :\n",
       "      original        bias    std. error\n",
       "t1* 39.9358610  0.0544513229 0.841289790\n",
       "t2* -0.1578447 -0.0006170901 0.007343073"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "boot(Auto, boot.fn, 1000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 70,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<caption>A matrix: 2 × 4 of type dbl</caption>\n",
       "<thead>\n",
       "\t<tr><th></th><th scope=col>Estimate</th><th scope=col>Std. Error</th><th scope=col>t value</th><th scope=col>Pr(&gt;|t|)</th></tr>\n",
       "</thead>\n",
       "<tbody>\n",
       "\t<tr><th scope=row>(Intercept)</th><td>39.9358610</td><td>0.717498656</td><td> 55.65984</td><td>1.220362e-187</td></tr>\n",
       "\t<tr><th scope=row>horsepower</th><td>-0.1578447</td><td>0.006445501</td><td>-24.48914</td><td> 7.031989e-81</td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "A matrix: 2 × 4 of type dbl\n",
       "\\begin{tabular}{r|llll}\n",
       "  & Estimate & Std. Error & t value & Pr(>\\textbar{}t\\textbar{})\\\\\n",
       "\\hline\n",
       "\t(Intercept) & 39.9358610 & 0.717498656 &  55.65984 & 1.220362e-187\\\\\n",
       "\thorsepower & -0.1578447 & 0.006445501 & -24.48914 &  7.031989e-81\\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "A matrix: 2 × 4 of type dbl\n",
       "\n",
       "| <!--/--> | Estimate | Std. Error | t value | Pr(&gt;|t|) |\n",
       "|---|---|---|---|---|\n",
       "| (Intercept) | 39.9358610 | 0.717498656 |  55.65984 | 1.220362e-187 |\n",
       "| horsepower | -0.1578447 | 0.006445501 | -24.48914 |  7.031989e-81 |\n",
       "\n"
      ],
      "text/plain": [
       "            Estimate   Std. Error  t value   Pr(>|t|)     \n",
       "(Intercept) 39.9358610 0.717498656  55.65984 1.220362e-187\n",
       "horsepower  -0.1578447 0.006445501 -24.48914  7.031989e-81"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "summary(lm(mpg~horsepower, data=Auto))$coef"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 72,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "\n",
       "ORDINARY NONPARAMETRIC BOOTSTRAP\n",
       "\n",
       "\n",
       "Call:\n",
       "boot(data = Auto, statistic = boot.fn, R = 1000)\n",
       "\n",
       "\n",
       "Bootstrap Statistics :\n",
       "        original        bias     std. error\n",
       "t1* 56.900099702  3.511640e-02 2.0300222526\n",
       "t2* -0.466189630 -7.080834e-04 0.0324241984\n",
       "t3*  0.001230536  2.840324e-06 0.0001172164"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "boot.fn = function(data, index) {\n",
    "    coefficients(lm(mpg~horsepower+I(horsepower^2), data=data, subset=index))\n",
    "}\n",
    "\n",
    "set.seed(1)\n",
    "boot(Auto, boot.fn, 1000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "R",
   "language": "R",
   "name": "ir"
  },
  "language_info": {
   "codemirror_mode": "r",
   "file_extension": ".r",
   "mimetype": "text/x-r-source",
   "name": "R",
   "pygments_lexer": "r",
   "version": "4.0.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}