lab_trees.ipynb

{
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3", 
      "name": "python3", 
      "language": "python"
    }, 
    "language_info": {
      "mimetype": "text/x-python", 
      "nbconvert_exporter": "python", 
      "name": "python", 
      "file_extension": ".py", 
      "version": "3.6.8", 
      "pygments_lexer": "ipython3", 
      "codemirror_mode": {
        "version": 3, 
        "name": "ipython"
      }
    }
  }, 
  "nbformat": 4, 
  "nbformat_minor": 2, 
  "cells": [
    {
      "source": [
        "# Trees\n", 
        "\n", 
        "(based on cs109a notebook by Rahul Dave, Kevin Raeder, and Pavlos Protopapas)"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "source": [
        "We will look here into the practicalities of fitting regression trees, random forests, and boosted trees. These involve out-of-bound estmates and cross-validation, and how you might want to deal with hyperparameters in these models. Along the way we will play a little bit with different loss functions, so that you start thinking about what goes in general into cooking up a machine learning model."
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": 21, 
      "cell_type": "code", 
      "source": [
        "%matplotlib inline\n", 
        "import numpy as np\n", 
        "import scipy as sp\n", 
        "import matplotlib as mpl\n", 
        "import matplotlib.cm as cm\n", 
        "import matplotlib.pyplot as plt\n", 
        "import pandas as pd\n", 
        "pd.set_option('display.width', 500)\n", 
        "pd.set_option('display.max_columns', 100)\n", 
        "pd.set_option('display.notebook_repr_html', True)\n", 
        "import seaborn.apionly as sns\n", 
        "sns.set_style(\"whitegrid\")\n", 
        "sns.set_context(\"poster\")"
      ], 
      "outputs": [], 
      "metadata": {}
    }, 
    {
      "source": [
        "## Dataset\n", 
        "\n", 
        "First, the data. We will be attempting to predict the presidential election results (at the county level) from 2016, measured as 'votergap' = (trump - clinton) in percentage points, based mostly on demographic features of those counties.  Let's quick take a peak at the data:"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": 22, 
      "cell_type": "code", 
      "source": [
        "elect_df = pd.read_csv(\"data/county_level_election.csv\")\n", 
        "elect_df.head()"
      ], 
      "outputs": [], 
      "metadata": {}
    }, 
    {
      "execution_count": 23, 
      "cell_type": "code", 
      "source": [
        "from sklearn.model_selection import train_test_split"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": true
      }
    }, 
    {
      "execution_count": 24, 
      "cell_type": "code", 
      "source": [
        "# split 80/20 train-test\n", 
        "Xdf = elect_df[['population','hispanic','minority','female','unemployed','income','nodegree','bachelor','inactivity','obesity','density','cancer']]\n", 
        "response = elect_df['votergap']\n", 
        "Xtraindf, Xtestdf, ytrain, ytest = train_test_split(Xdf, response, test_size=0.2, random_state=1983)\n"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": true
      }
    }, 
    {
      "execution_count": 25, 
      "cell_type": "code", 
      "source": [
        "plt.hist(ytrain)\n", 
        "Xtraindf.hist(column=['minority', 'population','hispanic','female']);"
      ], 
      "outputs": [], 
      "metadata": {}
    }, 
    {
      "source": [
        "## General Trees\n", 
        "\n", 
        "We could use a simple Decision Tree regressor to predict votergap. \n", 
        "\n", 
        "This is what you ought to keep in mind about decision trees.\n", 
        "\n", 
        "from the docs:\n", 
        "```\n", 
        "max_depth : int or None, optional (default=None)\n", 
        "The maximum depth of the tree. If None, then nodes are expanded until all leaves are pure or until all leaves contain less than min_samples_split samples.\n", 
        "min_samples_split : int, float, optional (default=2)\n", 
        "```\n", 
        "\n", 
        "- The deeper the tree, the more prone you are to overfitting.\n", 
        "- The smaller `min_samples_split`, the more the overfitting. \n", 
        "- One may use `min_samples_leaf` instead. More samples per leaf, the higher the bias."
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": 26, 
      "cell_type": "code", 
      "source": [
        "from sklearn.tree import DecisionTreeRegressor\n", 
        "#x = np.arange(0, 2*np.pi, 0.1)\n", 
        "#y = np.sin(x) + 0.1*np.random.normal(size=x.shape[0])\n", 
        "x = Xtraindf['minority'].values\n", 
        "o = np.argsort(x)\n", 
        "x = x[o]\n", 
        "y = ytrain.values\n", 
        "y = y[o]\n"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": true
      }
    }, 
    {
      "execution_count": 27, 
      "cell_type": "code", 
      "source": [
        "plt.plot(x,y, '.');"
      ], 
      "outputs": [], 
      "metadata": {}
    }, 
    {
      "execution_count": 28, 
      "cell_type": "code", 
      "source": [
        "plt.plot(np.log(x),y, '.');"
      ], 
      "outputs": [], 
      "metadata": {}
    }, 
    {
      "source": [
        "Which of the two versions of 'minority' would be a better choice to use as a predictor for inference?  For prediction?"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": 29, 
      "cell_type": "code", 
      "source": [
        "plt.plot(np.log(x),y,'.')\n", 
        "xx = np.log(x).reshape(-1,1)\n", 
        "for i in [1,2]:\n", 
        "    dtree = DecisionTreeRegressor(max_depth=i)\n", 
        "    dtree.fit(xx, y)\n", 
        "    plt.plot(np.log(x), dtree.predict(xx), label=str(i), alpha=1-i/10, lw=4)\n", 
        "plt.legend();"
      ], 
      "outputs": [], 
      "metadata": {}
    }, 
    {
      "execution_count": 30, 
      "cell_type": "code", 
      "source": [
        "plt.plot(np.log(x),y,'.')\n", 
        "xx = np.log(x).reshape(-1,1)\n", 
        "for i in [500,200,100,20]:\n", 
        "    dtree = DecisionTreeRegressor(min_samples_split=i)\n", 
        "    dtree.fit(xx, y)\n", 
        "    plt.plot(np.log(x), dtree.predict(xx), label=str(i), alpha=0.8, lw=4)\n", 
        "plt.legend();"
      ], 
      "outputs": [], 
      "metadata": {}
    }, 
    {
      "execution_count": 31, 
      "cell_type": "code", 
      "source": [
        "plt.plot(np.log(x),y,'.')\n", 
        "xx = np.log(x).reshape(-1,1)\n", 
        "for i in [500,200,100,20]:\n", 
        "    dtree = DecisionTreeRegressor(max_depth=6, min_samples_split=i)\n", 
        "    dtree.fit(xx, y)\n", 
        "    plt.plot(np.log(x), dtree.predict(xx), label=str(i), alpha=0.8, lw=4)\n", 
        "plt.legend();"
      ], 
      "outputs": [], 
      "metadata": {}
    }, 
    {
      "execution_count": 32, 
      "cell_type": "code", 
      "source": [
        "#let's also include logminority as a predictor going forward\n", 
        "xtemp = np.log(Xtraindf['minority'].values)\n", 
        "Xtraindf = Xtraindf.assign(logminority = xtemp)\n", 
        "Xtestdf = Xtestdf.assign(logminority = np.log(Xtestdf['minority'].values))\n", 
        "Xtraindf.head()\n"
      ], 
      "outputs": [], 
      "metadata": {}
    }, 
    {
      "source": [
        "# Exercise 1\n", 
        "\n", 
        "1. Perform 5-fold cross-validation to determine what the best `max_depth` would be for a single regression tree using the entire 'Xtrain' feature set. \n", 
        "2. Visualize the results with mean +/- 2 sd's across the validation sets."
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": 48, 
      "cell_type": "code", 
      "source": [
        "from sklearn.model_selection import cross_val_score\n", 
        "\n", 
        "depths = list(range(1, 21))\n", 
        "# train_scores = []\n", 
        "# cvmeans = []\n", 
        "# cvstds = []\n", 
        "# cv_scores = []\n", 
        "# your code here\n"
      ], 
      "outputs": [], 
      "metadata": {
        "collapsed": true
      }
    }, 
    {
      "execution_count": 49, 
      "cell_type": "code", 
      "source": [
        "# your code here\n"
      ], 
      "outputs": [], 
      "metadata": {}
    }, 
    {
      "source": [
        "Based on the plot above, it's clear that the training $R^2$ increases towards zero as `max_depth` increases, while the test set's $R^2$ maxes out around 62-63% around a `max_depth` of 7-9. Any of those values would be reasonable to choose for a best predictive model.  To prevent overfitting and for parsimony, Iwould choose a `max_depth` of 7.  "
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": 54, 
      "cell_type": "code", 
      "source": [
        "dtree = DecisionTreeRegressor(max_depth=6)\n", 
        "print(dtree.fit(Xtraindf, ytrain).score(Xtraindf, ytrain))\n", 
        "dtree.score(Xtestdf, ytest)"
      ], 
      "outputs": [], 
      "metadata": {}
    }, 
    {
      "source": [
        "Still some overfitting though"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "source": [
        "Trees are very nice models in that they give us feature importances:"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": 55, 
      "cell_type": "code", 
      "source": [
        "pd.Series(dtree.feature_importances_,index=list(Xtraindf)).sort_values().plot(kind=\"barh\")"
      ], 
      "outputs": [], 
      "metadata": {}
    }, 
    {
      "source": [
        "Ok with this discussion in mind, lets improve this model by Bagging."
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "source": [
        "## Bootstrap-Aggregating (called Bagging)\n", 
        "\n", 
        "Whats the basic idea?\n", 
        "\n", 
        "- A Single Decision tree is likely to overfit.\n", 
        "- So lets introduce replication through Bootstrap sampling.\n", 
        "- **Bagging** uses bootstrap resampling to create different training datasets. This way each training will give us a different tree.\n", 
        "- Added bonus: the left off points can be used to as a natural \"validation\" set, so no need to leave more data out for cross-validation\n", 
        "- Since we have many trees that we will **average over for prediction**, we can choose a large `max_depth` and we are ok as we will rely on the law of large numbers to shrink this large variance, low bias approach for each individual tree."
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": 36, 
      "cell_type": "code", 
      "source": [
        "from sklearn.utils import resample\n", 
        "\n", 
        "ntrees = 500\n", 
        "estimators = []\n", 
        "R2s = []\n", 
        "yhats_test = np.zeros((Xtestdf.shape[0], ntrees))\n", 
        "\n", 
        "plt.plot(np.log(x),y,'.')\n", 
        "for i in range(ntrees):\n", 
        "    simpletree = DecisionTreeRegressor(max_depth=3)\n", 
        "    boot_xx, boot_y = resample(Xtraindf[['logminority']], ytrain)\n", 
        "    estimators = np.append(estimators,simpletree.fit(boot_xx, boot_y))\n", 
        "    R2s = np.append(R2s,simpletree.score(Xtestdf[['logminority']], ytest))\n", 
        "    yhats_test[:,i] = simpletree.predict(Xtestdf[['logminority']])\n", 
        "    plt.plot(np.log(x), simpletree.predict(np.log(x).reshape(-1,1)), 'red', alpha=0.05)"
      ], 
      "outputs": [], 
      "metadata": {}
    }, 
    {
      "execution_count": 37, 
      "cell_type": "code", 
      "source": [
        "yhats_test.shape"
      ], 
      "outputs": [], 
      "metadata": {}
    }, 
    {
      "source": [
        "# Exercise\n", 
        "\n", 
        "1. Edit the code below (which is just copied from above) to refit many bagged trees on the entire `Xtrain` feature set (without the plot...lots of predictors now so difficult to plot). \n", 
        "2. Summarize how each of the separate trees performed (both numerically and visually) using $R^2$ as the metric.  How do they perform on average?\n", 
        "3. Combine the trees into one prediction and evaluate it using $R^2$.\n", 
        "4. Briefly discuss the results.  How will the results above change if 'max_depth' is increased?  What if it is decreased?"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": 38, 
      "cell_type": "code", 
      "source": [
        "from sklearn.metrics import r2_score\n", 
        "\n", 
        "\n", 
        "ntrees = 500\n", 
        "estimators = []\n", 
        "R2s = []\n", 
        "yhats_test = np.zeros((Xtestdf.shape[0], ntrees))\n", 
        "\n", 
        "# your code here\n"
      ], 
      "outputs": [], 
      "metadata": {}
    }, 
    {
      "source": [
        "4. Notice that the aggregated results, by averaging the predictions across the 500 different Regression Tree models, scores and $R^2$ in the test set of 64.9% which outperforms all of the separate individual models.  We would expect this to improve further for larger values of `max_depth`, and in fact should imrpove quite a lot.  This Bagging improves the variability of 'high variance' models by leveraging the law of large numbers, but does not improve bias nearly as much.  Essentially, it is an indirect way of 'smoothing' these discretized step function by essentially jittering where those jumps occur.  See the visual above.  \n", 
        "\n", 
        "Try a max_depth of 15"
      ], 
      "cell_type": "markdown", 
      "metadata": {}
    }, 
    {
      "execution_count": 39, 
      "cell_type": "code", 
      "source": [
        "from sklearn.metrics import r2_score\n", 
        "\n", 
        "\n", 
        "# your code here\n"
      ], 
      "outputs": [], 
      "metadata": {}
    }
  ]
}