nb2_sklearn-lr-rf.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Notebook 2: Linear regression and random forests with scikit-learn\n",
    "\n",
    "Now that we have prepared our data we can start building machine learning models.\n",
    "\n",
    "We will start with the simplest model: Linear regression"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import pickle\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import seaborn as sns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "DATADIR = './dataset1/'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Copy from previous notebook\n",
    "def create_sub(preds, fn=None):\n",
    "    df =  pd.DataFrame({'Id': range(len(preds)), 'Expected': preds})\n",
    "    if fn is not None: df.to_csv(DATADIR + fn, index=False)\n",
    "    return df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Load the preprocessed data\n",
    "with open('./tmp/preproc_data.pkl', 'rb') as f:\n",
    "    X_train, y_train, X_valid, y_valid, X_test = pickle.load(f)\n",
    "with open('./tmp/dfs.pkl', 'rb') as f:\n",
    "    df_train, df_test = pickle.load(f)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Linear regression\n",
    "\n",
    "$y = \\sum_i a_i x_i + b$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Compute baseline mse\n",
    "\n",
    "Before we start let's check which score we have to beat."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def mse(y_true, y_pred): return ((y_true - y_pred)**2).mean()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "4.135858792431942"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mse(y_valid, df_train.t2m_fc_mean[df_train.time >= '2015-01-01'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### sklearn API\n",
    "\n",
    "All sklearn models work in exactly the same way."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from sklearn.linear_model import LinearRegression"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Create the model object\n",
    "lr = LinearRegression()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "LinearRegression(copy_X=True, fit_intercept=True, n_jobs=1, normalize=False)"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Fit the model parameters with the training data\n",
    "lr.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now the parameters $a_i$ and $b$ have been fit to our training data. \n",
    "\n",
    "We can not look at the individual parameters."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Index(['t2m_fc_mean', 'orog', 'station_alt', 'station_lat', 'station_lon',\n",
       "       'u_pl500_fc_mean', 'v_pl500_fc_mean', 'gh_pl500_fc_mean',\n",
       "       'u_pl850_fc_mean', 'v_pl850_fc_mean', 'q_pl850_fc_mean', 'cape_fc_mean',\n",
       "       'sp_fc_mean', 'tcc_fc_mean', 'sshf_fc_mean', 'slhf_fc_mean',\n",
       "       'u10_fc_mean', 'v10_fc_mean', 'ssr_fc_mean', 'str_fc_mean',\n",
       "       'd2m_fc_mean', 'sm_fc_mean'],\n",
       "      dtype='object')"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X_train.columns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 5.20734246,  0.53522914, -0.9193415 , -0.14865423, -0.10651631,\n",
       "        0.01482758, -0.06200976,  0.22355911,  0.08053755, -0.01852129,\n",
       "        0.14741519, -0.03510716, -0.17207246,  0.17164729,  0.54843072,\n",
       "       -0.03918102, -0.06724854,  0.13348599,  0.77616008,  0.21138807,\n",
       "        0.61088598, -0.04088925])"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lr.coef_"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "7.448795068185583"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lr.intercept_"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Wouldn't you expect the coefficient for `t2m_fc_mean` to be around 0 and the bias to be close to one? Why is this not the case?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Make predictions for validation set\n",
    "preds = lr.predict(X_valid)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3.2422963614329663"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mse(preds, y_valid)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Id</th>\n",
       "      <th>Expected</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0</td>\n",
       "      <td>4.195768</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>1</td>\n",
       "      <td>1.528290</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2</td>\n",
       "      <td>1.475483</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>3</td>\n",
       "      <td>4.318616</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>4</td>\n",
       "      <td>2.337542</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   Id  Expected\n",
       "0   0  4.195768\n",
       "1   1  1.528290\n",
       "2   2  1.475483\n",
       "3   3  4.318616\n",
       "4   4  2.337542"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Submit to Kaggle\n",
    "df_sub = create_sub(lr.predict(X_test), 'lr.csv'); df_sub.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### The R$^2$ score\n",
    "\n",
    "The R$^2$ score is a common skill metric in machine learning: https://en.wikipedia.org/wiki/Coefficient_of_determination"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.9187372901848857"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lr.score(X_valid, y_valid)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def print_scores(m, X_train=X_train, X_valid=X_valid):\n",
    "    print('Train R2 = ', m.score(X_train, y_train), \n",
    "          ', Valid R2 = ', m.score(X_valid, y_valid), ', Valid MSE = ', \n",
    "          mse(m.predict(X_valid), y_valid))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Train R2 =  0.931943086657 , Valid R2 =  0.918737290185 , Valid MSE =  3.2422963614329663\n"
     ]
    }
   ],
   "source": [
    "print_scores(lr)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can see that the training and validation score are roughly similar, indicating that the linear regression doesn't seem to be very prone to overfitting."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Decicion trees\n",
    "\n",
    "Linear regression improves the score but we can do better. One very common machine learning algorithm are decision trees. They are capable of representing non-linear relationships. In scikit-learn we can use the exact same API as before."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.tree import DecisionTreeRegressor\n",
    "from sklearn.tree import export_graphviz\n",
    "from IPython.display import display\n",
    "import IPython, graphviz, re"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualize a tree"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def draw_tree(t, df, size=10, ratio=0.6, precision=0):\n",
    "    \"\"\" From fast.ai machine learning course.\n",
    "    Draws a representation of a random forest in IPython.\n",
    "    Parameters:\n",
    "    -----------\n",
    "    t: The tree you wish to draw\n",
    "    df: The data used to train the tree. This is used to get the names of the features.\n",
    "    \"\"\"\n",
    "    s=export_graphviz(t, out_file=None, feature_names=df.columns, filled=True,\n",
    "                      special_characters=True, rotate=True, precision=precision)\n",
    "    IPython.display.display(graphviz.Source(re.sub('Tree {',\n",
    "       f'Tree {{ size={size}; ratio={ratio}', s)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's start with a very shallow tree which we can visualize. Again we are using the exact same sklearn API."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "dt = DecisionTreeRegressor(max_depth=3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Wall time: 5.33 s\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "DecisionTreeRegressor(criterion='mse', max_depth=3, max_features=None,\n",
       "           max_leaf_nodes=None, min_impurity_decrease=0.0,\n",
       "           min_impurity_split=None, min_samples_leaf=1,\n",
       "           min_samples_split=2, min_weight_fraction_leaf=0.0,\n",
       "           presort=False, random_state=None, splitter='best')"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "%time dt.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/svg+xml": [
       "<?xml version=\"1.0\" encoding=\"UTF-8\" standalone=\"no\"?>\r\n",
       "<!DOCTYPE svg PUBLIC \"-//W3C//DTD SVG 1.1//EN\"\r\n",
       " \"http://www.w3.org/Graphics/SVG/1.1/DTD/svg11.dtd\">\r\n",
       "<!-- Generated by graphviz version 2.38.0 (20140413.2041)\r\n",
       " -->\r\n",
       "<!-- Title: Tree Pages: 1 -->\r\n",
       "<svg width=\"720pt\" height=\"434pt\"\r\n",
       " viewBox=\"0.00 0.00 720.00 434.49\" xmlns=\"http://www.w3.org/2000/svg\" xmlns:xlink=\"http://www.w3.org/1999/xlink\">\r\n",
       "<g id=\"graph0\" class=\"graph\" transform=\"scale(0.778659 0.778659) rotate(0) translate(4 554)\">\r\n",
       "<title>Tree</title>\r\n",
       "<polygon fill=\"white\" stroke=\"none\" points=\"-4,4 -4,-554 920.667,-554 920.667,4 -4,4\"/>\r\n",
       "<!-- 0 -->\r\n",
       "<g id=\"node1\" class=\"node\"><title>0</title>\r\n",
       "<polygon fill=\"#e58139\" fill-opacity=\"0.588235\" stroke=\"black\" points=\"161.667,-308.5 31.6667,-308.5 31.6667,-240.5 161.667,-240.5 161.667,-308.5\"/>\r\n",
       "<text text-anchor=\"start\" x=\"39.6667\" y=\"-293.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">t2m_fc_mean ≤ &#45;0.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"65.1667\" y=\"-278.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">mse = 44.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"45.6667\" y=\"-263.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">samples = 728008</text>\r\n",
       "<text text-anchor=\"start\" x=\"65.6667\" y=\"-248.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">value = 7.0</text>\r\n",
       "</g>\r\n",
       "<!-- 1 -->\r\n",
       "<g id=\"node2\" class=\"node\"><title>1</title>\r\n",
       "<polygon fill=\"#e58139\" fill-opacity=\"0.356863\" stroke=\"black\" points=\"408.667,-361.5 278.667,-361.5 278.667,-293.5 408.667,-293.5 408.667,-361.5\"/>\r\n",
       "<text text-anchor=\"start\" x=\"286.667\" y=\"-346.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">t2m_fc_mean ≤ &#45;1.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"312.167\" y=\"-331.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">mse = 19.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"292.667\" y=\"-316.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">samples = 333362</text>\r\n",
       "<text text-anchor=\"start\" x=\"312.667\" y=\"-301.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">value = 2.0</text>\r\n",
       "</g>\r\n",
       "<!-- 0&#45;&gt;1 -->\r\n",
       "<g id=\"edge1\" class=\"edge\"><title>0&#45;&gt;1</title>\r\n",
       "<path fill=\"none\" stroke=\"black\" d=\"M161.718,-288.356C194.494,-295.446 234.563,-304.114 268.656,-311.489\"/>\r\n",
       "<polygon fill=\"black\" stroke=\"black\" points=\"268.051,-314.939 278.565,-313.633 269.531,-308.098 268.051,-314.939\"/>\r\n",
       "<text text-anchor=\"middle\" x=\"257.549\" y=\"-323.473\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">True</text>\r\n",
       "</g>\r\n",
       "<!-- 8 -->\r\n",
       "<g id=\"node9\" class=\"node\"><title>8</title>\r\n",
       "<polygon fill=\"#e58139\" fill-opacity=\"0.780392\" stroke=\"black\" points=\"406.167,-255.5 281.167,-255.5 281.167,-187.5 406.167,-187.5 406.167,-255.5\"/>\r\n",
       "<text text-anchor=\"start\" x=\"289.167\" y=\"-240.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">t2m_fc_mean ≤ 1.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"312.167\" y=\"-225.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">mse = 14.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"292.667\" y=\"-210.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">samples = 394646</text>\r\n",
       "<text text-anchor=\"start\" x=\"309.167\" y=\"-195.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">value = 12.0</text>\r\n",
       "</g>\r\n",
       "<!-- 0&#45;&gt;8 -->\r\n",
       "<g id=\"edge8\" class=\"edge\"><title>0&#45;&gt;8</title>\r\n",
       "<path fill=\"none\" stroke=\"black\" d=\"M161.718,-260.644C195.268,-253.386 236.46,-244.475 271.062,-236.99\"/>\r\n",
       "<polygon fill=\"black\" stroke=\"black\" points=\"272.068,-240.353 281.102,-234.818 270.588,-233.512 272.068,-240.353\"/>\r\n",
       "<text text-anchor=\"middle\" x=\"260.087\" y=\"-217.578\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">False</text>\r\n",
       "</g>\r\n",
       "<!-- 2 -->\r\n",
       "<g id=\"node3\" class=\"node\"><title>2</title>\r\n",
       "<polygon fill=\"#e58139\" fill-opacity=\"0.176471\" stroke=\"black\" points=\"655.667,-486.5 525.667,-486.5 525.667,-418.5 655.667,-418.5 655.667,-486.5\"/>\r\n",
       "<text text-anchor=\"start\" x=\"533.667\" y=\"-471.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">t2m_fc_mean ≤ &#45;2.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"559.167\" y=\"-456.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">mse = 15.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"539.667\" y=\"-441.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">samples = 108193</text>\r\n",
       "<text text-anchor=\"start\" x=\"557.167\" y=\"-426.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">value = &#45;3.0</text>\r\n",
       "</g>\r\n",
       "<!-- 1&#45;&gt;2 -->\r\n",
       "<g id=\"edge2\" class=\"edge\"><title>1&#45;&gt;2</title>\r\n",
       "<path fill=\"none\" stroke=\"black\" d=\"M408.718,-360.179C441.775,-377.045 482.251,-397.696 516.532,-415.187\"/>\r\n",
       "<polygon fill=\"black\" stroke=\"black\" points=\"515.066,-418.368 525.565,-419.795 518.248,-412.132 515.066,-418.368\"/>\r\n",
       "</g>\r\n",
       "<!-- 5 -->\r\n",
       "<g id=\"node6\" class=\"node\"><title>5</title>\r\n",
       "<polygon fill=\"#e58139\" fill-opacity=\"0.447059\" stroke=\"black\" points=\"655.667,-361.5 525.667,-361.5 525.667,-293.5 655.667,-293.5 655.667,-361.5\"/>\r\n",
       "<text text-anchor=\"start\" x=\"533.667\" y=\"-346.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">t2m_fc_mean ≤ &#45;1.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"562.667\" y=\"-331.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">mse = 7.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"539.667\" y=\"-316.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">samples = 225169</text>\r\n",
       "<text text-anchor=\"start\" x=\"559.667\" y=\"-301.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">value = 4.0</text>\r\n",
       "</g>\r\n",
       "<!-- 1&#45;&gt;5 -->\r\n",
       "<g id=\"edge5\" class=\"edge\"><title>1&#45;&gt;5</title>\r\n",
       "<path fill=\"none\" stroke=\"black\" d=\"M408.718,-327.5C441.353,-327.5 481.219,-327.5 515.216,-327.5\"/>\r\n",
       "<polygon fill=\"black\" stroke=\"black\" points=\"515.565,-331 525.565,-327.5 515.565,-324 515.565,-331\"/>\r\n",
       "</g>\r\n",
       "<!-- 3 -->\r\n",
       "<g id=\"node4\" class=\"node\"><title>3</title>\r\n",
       "<polygon fill=\"none\" stroke=\"black\" points=\"884.167,-550 773.167,-550 773.167,-497 884.167,-497 884.167,-550\"/>\r\n",
       "<text text-anchor=\"start\" x=\"797.167\" y=\"-534.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">mse = 17.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"781.167\" y=\"-519.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">samples = 27937</text>\r\n",
       "<text text-anchor=\"start\" x=\"795.167\" y=\"-504.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">value = &#45;7.0</text>\r\n",
       "</g>\r\n",
       "<!-- 2&#45;&gt;3 -->\r\n",
       "<g id=\"edge3\" class=\"edge\"><title>2&#45;&gt;3</title>\r\n",
       "<path fill=\"none\" stroke=\"black\" d=\"M655.843,-471.807C689.261,-481.861 730.008,-494.12 763.434,-504.176\"/>\r\n",
       "<polygon fill=\"black\" stroke=\"black\" points=\"762.521,-507.556 773.105,-507.085 764.537,-500.853 762.521,-507.556\"/>\r\n",
       "</g>\r\n",
       "<!-- 4 -->\r\n",
       "<g id=\"node5\" class=\"node\"><title>4</title>\r\n",
       "<polygon fill=\"#e58139\" fill-opacity=\"0.239216\" stroke=\"black\" points=\"884.167,-479 773.167,-479 773.167,-426 884.167,-426 884.167,-479\"/>\r\n",
       "<text text-anchor=\"start\" x=\"800.667\" y=\"-463.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">mse = 5.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"781.167\" y=\"-448.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">samples = 80256</text>\r\n",
       "<text text-anchor=\"start\" x=\"795.167\" y=\"-433.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">value = &#45;1.0</text>\r\n",
       "</g>\r\n",
       "<!-- 2&#45;&gt;4 -->\r\n",
       "<g id=\"edge4\" class=\"edge\"><title>2&#45;&gt;4</title>\r\n",
       "<path fill=\"none\" stroke=\"black\" d=\"M655.843,-452.5C689.117,-452.5 729.659,-452.5 763.004,-452.5\"/>\r\n",
       "<polygon fill=\"black\" stroke=\"black\" points=\"763.105,-456 773.105,-452.5 763.105,-449 763.105,-456\"/>\r\n",
       "</g>\r\n",
       "<!-- 6 -->\r\n",
       "<g id=\"node7\" class=\"node\"><title>6</title>\r\n",
       "<polygon fill=\"#e58139\" fill-opacity=\"0.384314\" stroke=\"black\" points=\"887.667,-408 769.667,-408 769.667,-355 887.667,-355 887.667,-408\"/>\r\n",
       "<text text-anchor=\"start\" x=\"800.667\" y=\"-392.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">mse = 4.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"777.667\" y=\"-377.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">samples = 112684</text>\r\n",
       "<text text-anchor=\"start\" x=\"797.667\" y=\"-362.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">value = 2.0</text>\r\n",
       "</g>\r\n",
       "<!-- 5&#45;&gt;6 -->\r\n",
       "<g id=\"edge6\" class=\"edge\"><title>5&#45;&gt;6</title>\r\n",
       "<path fill=\"none\" stroke=\"black\" d=\"M655.843,-342.184C688.004,-349.543 726.954,-358.456 759.633,-365.933\"/>\r\n",
       "<polygon fill=\"black\" stroke=\"black\" points=\"759.03,-369.385 769.559,-368.204 760.592,-362.562 759.03,-369.385\"/>\r\n",
       "</g>\r\n",
       "<!-- 7 -->\r\n",
       "<g id=\"node8\" class=\"node\"><title>7</title>\r\n",
       "<polygon fill=\"#e58139\" fill-opacity=\"0.509804\" stroke=\"black\" points=\"887.667,-337 769.667,-337 769.667,-284 887.667,-284 887.667,-337\"/>\r\n",
       "<text text-anchor=\"start\" x=\"800.667\" y=\"-321.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">mse = 4.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"777.667\" y=\"-306.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">samples = 112485</text>\r\n",
       "<text text-anchor=\"start\" x=\"797.667\" y=\"-291.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">value = 5.0</text>\r\n",
       "</g>\r\n",
       "<!-- 5&#45;&gt;7 -->\r\n",
       "<g id=\"edge7\" class=\"edge\"><title>5&#45;&gt;7</title>\r\n",
       "<path fill=\"none\" stroke=\"black\" d=\"M655.843,-322.877C687.865,-320.57 726.618,-317.779 759.21,-315.431\"/>\r\n",
       "<polygon fill=\"black\" stroke=\"black\" points=\"759.836,-318.895 769.559,-314.686 759.333,-311.913 759.836,-318.895\"/>\r\n",
       "</g>\r\n",
       "<!-- 9 -->\r\n",
       "<g id=\"node10\" class=\"node\"><title>9</title>\r\n",
       "<polygon fill=\"#e58139\" fill-opacity=\"0.682353\" stroke=\"black\" points=\"653.167,-255.5 528.167,-255.5 528.167,-187.5 653.167,-187.5 653.167,-255.5\"/>\r\n",
       "<text text-anchor=\"start\" x=\"536.167\" y=\"-240.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">t2m_fc_mean ≤ 0.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"562.667\" y=\"-225.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">mse = 6.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"539.667\" y=\"-210.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">samples = 218471</text>\r\n",
       "<text text-anchor=\"start\" x=\"556.167\" y=\"-195.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">value = 10.0</text>\r\n",
       "</g>\r\n",
       "<!-- 8&#45;&gt;9 -->\r\n",
       "<g id=\"edge9\" class=\"edge\"><title>8&#45;&gt;9</title>\r\n",
       "<path fill=\"none\" stroke=\"black\" d=\"M406.174,-221.5C440.104,-221.5 482.415,-221.5 517.854,-221.5\"/>\r\n",
       "<polygon fill=\"black\" stroke=\"black\" points=\"518.134,-225 528.134,-221.5 518.134,-218 518.134,-225\"/>\r\n",
       "</g>\r\n",
       "<!-- 12 -->\r\n",
       "<g id=\"node13\" class=\"node\"><title>12</title>\r\n",
       "<polygon fill=\"#e58139\" fill-opacity=\"0.898039\" stroke=\"black\" points=\"653.167,-131.5 528.167,-131.5 528.167,-63.5 653.167,-63.5 653.167,-131.5\"/>\r\n",
       "<text text-anchor=\"start\" x=\"536.167\" y=\"-116.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">t2m_fc_mean ≤ 1.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"562.667\" y=\"-101.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">mse = 7.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"539.667\" y=\"-86.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">samples = 176175</text>\r\n",
       "<text text-anchor=\"start\" x=\"556.167\" y=\"-71.3\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">value = 15.0</text>\r\n",
       "</g>\r\n",
       "<!-- 8&#45;&gt;12 -->\r\n",
       "<g id=\"edge12\" class=\"edge\"><title>8&#45;&gt;12</title>\r\n",
       "<path fill=\"none\" stroke=\"black\" d=\"M406.174,-190.37C440.395,-173.05 483.142,-151.414 518.765,-133.385\"/>\r\n",
       "<polygon fill=\"black\" stroke=\"black\" points=\"520.792,-136.282 528.134,-128.643 517.631,-130.036 520.792,-136.282\"/>\r\n",
       "</g>\r\n",
       "<!-- 10 -->\r\n",
       "<g id=\"node11\" class=\"node\"><title>10</title>\r\n",
       "<polygon fill=\"#e58139\" fill-opacity=\"0.623529\" stroke=\"black\" points=\"887.667,-266 769.667,-266 769.667,-213 887.667,-213 887.667,-266\"/>\r\n",
       "<text text-anchor=\"start\" x=\"800.667\" y=\"-250.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">mse = 4.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"777.667\" y=\"-235.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">samples = 100215</text>\r\n",
       "<text text-anchor=\"start\" x=\"797.667\" y=\"-220.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">value = 8.0</text>\r\n",
       "</g>\r\n",
       "<!-- 9&#45;&gt;10 -->\r\n",
       "<g id=\"edge10\" class=\"edge\"><title>9&#45;&gt;10</title>\r\n",
       "<path fill=\"none\" stroke=\"black\" d=\"M653.365,-226.206C685.916,-228.688 725.933,-231.741 759.411,-234.294\"/>\r\n",
       "<polygon fill=\"black\" stroke=\"black\" points=\"759.338,-237.799 769.575,-235.069 759.87,-230.819 759.338,-237.799\"/>\r\n",
       "</g>\r\n",
       "<!-- 11 -->\r\n",
       "<g id=\"node12\" class=\"node\"><title>11</title>\r\n",
       "<polygon fill=\"#e58139\" fill-opacity=\"0.733333\" stroke=\"black\" points=\"887.667,-195 769.667,-195 769.667,-142 887.667,-142 887.667,-195\"/>\r\n",
       "<text text-anchor=\"start\" x=\"800.667\" y=\"-179.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">mse = 4.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"777.667\" y=\"-164.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">samples = 118256</text>\r\n",
       "<text text-anchor=\"start\" x=\"794.167\" y=\"-149.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">value = 11.0</text>\r\n",
       "</g>\r\n",
       "<!-- 9&#45;&gt;11 -->\r\n",
       "<g id=\"edge11\" class=\"edge\"><title>9&#45;&gt;11</title>\r\n",
       "<path fill=\"none\" stroke=\"black\" d=\"M653.365,-207.644C685.916,-200.334 725.933,-191.347 759.411,-183.829\"/>\r\n",
       "<polygon fill=\"black\" stroke=\"black\" points=\"760.585,-187.152 769.575,-181.546 759.051,-180.322 760.585,-187.152\"/>\r\n",
       "</g>\r\n",
       "<!-- 13 -->\r\n",
       "<g id=\"node14\" class=\"node\"><title>13</title>\r\n",
       "<polygon fill=\"#e58139\" fill-opacity=\"0.850980\" stroke=\"black\" points=\"887.667,-124 769.667,-124 769.667,-71 887.667,-71 887.667,-124\"/>\r\n",
       "<text text-anchor=\"start\" x=\"800.667\" y=\"-108.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">mse = 4.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"777.667\" y=\"-93.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">samples = 120810</text>\r\n",
       "<text text-anchor=\"start\" x=\"794.167\" y=\"-78.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">value = 14.0</text>\r\n",
       "</g>\r\n",
       "<!-- 12&#45;&gt;13 -->\r\n",
       "<g id=\"edge13\" class=\"edge\"><title>12&#45;&gt;13</title>\r\n",
       "<path fill=\"none\" stroke=\"black\" d=\"M653.365,-97.5C685.916,-97.5 725.933,-97.5 759.411,-97.5\"/>\r\n",
       "<polygon fill=\"black\" stroke=\"black\" points=\"759.575,-101 769.575,-97.5 759.575,-94.0001 759.575,-101\"/>\r\n",
       "</g>\r\n",
       "<!-- 14 -->\r\n",
       "<g id=\"node15\" class=\"node\"><title>14</title>\r\n",
       "<polygon fill=\"#e58139\" stroke=\"black\" points=\"884.167,-53 773.167,-53 773.167,-0 884.167,-0 884.167,-53\"/>\r\n",
       "<text text-anchor=\"start\" x=\"800.667\" y=\"-37.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">mse = 5.0</text>\r\n",
       "<text text-anchor=\"start\" x=\"781.167\" y=\"-22.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">samples = 55365</text>\r\n",
       "<text text-anchor=\"start\" x=\"794.167\" y=\"-7.8\" font-family=\"Times New Roman,serif\" font-size=\"14.00\">value = 18.0</text>\r\n",
       "</g>\r\n",
       "<!-- 12&#45;&gt;14 -->\r\n",
       "<g id=\"edge14\" class=\"edge\"><title>12&#45;&gt;14</title>\r\n",
       "<path fill=\"none\" stroke=\"black\" d=\"M653.365,-78.9383C687.175,-68.7667 729.039,-56.1718 763.265,-45.875\"/>\r\n",
       "<polygon fill=\"black\" stroke=\"black\" points=\"764.595,-49.1301 773.162,-42.8975 762.578,-42.4269 764.595,-49.1301\"/>\r\n",
       "</g>\r\n",
       "</g>\r\n",
       "</svg>\r\n"
      ],
      "text/plain": [
       "<graphviz.files.Source at 0x214315b0048>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "draw_tree(dt, df_train.drop(['t2m_obs', 'time', 'station'], 1))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here we can see that the splits are entirely based on `t2m_fc_mean` which is not surpising."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Build a full tree\n",
    "\n",
    "Now let's actually build a complicated tree."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "dt = DecisionTreeRegressor()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Wall time: 30.9 s\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "DecisionTreeRegressor(criterion='mse', max_depth=None, max_features=None,\n",
       "           max_leaf_nodes=None, min_impurity_decrease=0.0,\n",
       "           min_impurity_split=None, min_samples_leaf=1,\n",
       "           min_samples_split=2, min_weight_fraction_leaf=0.0,\n",
       "           presort=False, random_state=None, splitter='best')"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "%time dt.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Train R2 =  1.0 , Valid R2 =  0.843798011081 , Valid MSE =  6.232294510890336\n"
     ]
    }
   ],
   "source": [
    "print_scores(dt)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1291959"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "dt.tree_.node_count"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So we have a perfect score on the training dataset, but a pretty awful score on the validation set. This is an extreme case of overfitting!\n",
    "\n",
    "But this is to be expected. We built our tree to perfectly fit the training set, so that only one sample remains in the final leafs."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Reduce the complexity of the tree\n",
    "\n",
    "To combat this we can set the minimum number of samples in each final leaf to a higher value."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "dt = DecisionTreeRegressor(min_samples_leaf=200)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Wall time: 20.2 s\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "DecisionTreeRegressor(criterion='mse', max_depth=None, max_features=None,\n",
       "           max_leaf_nodes=None, min_impurity_decrease=0.0,\n",
       "           min_impurity_split=None, min_samples_leaf=200,\n",
       "           min_samples_split=2, min_weight_fraction_leaf=0.0,\n",
       "           presort=False, random_state=None, splitter='best')"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "%time dt.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Train R2 =  0.939008732906 , Valid R2 =  0.911272929509 , Valid MSE =  3.540116474929168\n"
     ]
    }
   ],
   "source": [
    "print_scores(dt)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The overfitting is less bad now but we still cannot get a better score than with linear regression. We need to do better!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Random forests\n",
    "\n",
    "Here is where random forests are coming in. As the name suggests, a random forest is an ensemble of decision trees, where each tree is build with a random subsample of the features and training data. This should help to prevent overfitting."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from sklearn.ensemble import RandomForestRegressor"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Fit a sample model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we are fitting several decision trees. This can be done in parallel by several processors. `n_jobs = -1` uses all available processors."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "rf = RandomForestRegressor(n_estimators=10, n_jobs=-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Wall time: 1min 27s\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "RandomForestRegressor(bootstrap=True, criterion='mse', max_depth=None,\n",
       "           max_features='auto', max_leaf_nodes=None,\n",
       "           min_impurity_decrease=0.0, min_impurity_split=None,\n",
       "           min_samples_leaf=1, min_samples_split=2,\n",
       "           min_weight_fraction_leaf=0.0, n_estimators=10, n_jobs=-1,\n",
       "           oob_score=False, random_state=None, verbose=0, warm_start=False)"
      ]
     },
     "execution_count": 34,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "%time rf.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Train R2 =  0.992415062421 , Valid R2 =  0.917427004519 , Valid MSE =  3.294575376695344\n"
     ]
    }
   ],
   "source": [
    "print_scores(rf)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Again, we see that we are overfitting quite badly. As before, we can tune the hyper-parameters to prevent overfitting"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Hyperparameters\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### min_samples_leaf\n",
    "\n",
    "The minimum number of sample is each leaf."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "rf = RandomForestRegressor(n_estimators=10, n_jobs=-1, min_samples_leaf=100)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Wall time: 1min 9s\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "RandomForestRegressor(bootstrap=True, criterion='mse', max_depth=None,\n",
       "           max_features='auto', max_leaf_nodes=None,\n",
       "           min_impurity_decrease=0.0, min_impurity_split=None,\n",
       "           min_samples_leaf=100, min_samples_split=2,\n",
       "           min_weight_fraction_leaf=0.0, n_estimators=10, n_jobs=-1,\n",
       "           oob_score=False, random_state=None, verbose=0, warm_start=False)"
      ]
     },
     "execution_count": 37,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "%time rf.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Train R2 =  0.946901912503 , Valid R2 =  0.919607573801 , Valid MSE =  3.2075729635844925\n"
     ]
    }
   ],
   "source": [
    "print_scores(rf)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### n_estimators\n",
    "\n",
    "The number of trees. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "rf = RandomForestRegressor(n_estimators=40, n_jobs=-1, min_samples_leaf=10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Wall time: 7min 16s\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "RandomForestRegressor(bootstrap=True, criterion='mse', max_depth=None,\n",
       "           max_features='auto', max_leaf_nodes=None,\n",
       "           min_impurity_decrease=0.0, min_impurity_split=None,\n",
       "           min_samples_leaf=10, min_samples_split=2,\n",
       "           min_weight_fraction_leaf=0.0, n_estimators=40, n_jobs=-1,\n",
       "           oob_score=False, random_state=None, verbose=0, warm_start=False)"
      ]
     },
     "execution_count": 40,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "%time rf.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Train R2 =  0.975899932589 , Valid R2 =  0.924019832822 , Valid MSE =  3.0315284851041664\n"
     ]
    }
   ],
   "source": [
    "print_scores(rf)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We are still getting a much better score on the training set. Partly we are still overfitting, but partly this is to be expected from our dataset and the way we chose the validation set. Can you explain why?\n",
    "\n",
    "Let's submit our new best model to kaggle."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Id</th>\n",
       "      <th>Expected</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0</td>\n",
       "      <td>3.934673</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>1</td>\n",
       "      <td>1.452402</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2</td>\n",
       "      <td>0.716910</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>3</td>\n",
       "      <td>3.903693</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>4</td>\n",
       "      <td>2.480547</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   Id  Expected\n",
       "0   0  3.934673\n",
       "1   1  1.452402\n",
       "2   2  0.716910\n",
       "3   3  3.903693\n",
       "4   4  2.480547"
      ]
     },
     "execution_count": 42,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Submit to Kaggle\n",
    "df_sub = create_sub(rf.predict(X_test), 'rf.csv'); df_sub.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Feature importance\n",
    "\n",
    "Getting a good prediction is nice, but ideally we would also understand why we are getting the prediction.\n",
    "\n",
    "For this we can look at feature importance. In feature importance, you randomly shuffle each column of the dataset and check by how much the prediction skill decreases."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 0.93886066,  0.00189192,  0.00730978,  0.00474286,  0.0022422 ,\n",
       "        0.00216329,  0.00212882,  0.00319009,  0.00231558,  0.00209869,\n",
       "        0.0032293 ,  0.00195781,  0.00137903,  0.00264125,  0.00224996,\n",
       "        0.00175749,  0.00234605,  0.00209426,  0.00188215,  0.0021487 ,\n",
       "        0.01009336,  0.00127674])"
      ]
     },
     "execution_count": 47,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "rf.feature_importances_"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Create a new dataframe for easy plotting\n",
    "fi_df = pd.DataFrame(\n",
    "    data={'Feature': X_train.columns, 'Feature importance': rf.feature_importances_}, \n",
    "    columns = ['Feature', 'Feature importance']\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "fi_df.sort_values('Feature importance', inplace=True, ascending=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 62,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\rasps\\Anaconda3\\lib\\site-packages\\seaborn\\categorical.py:1460: FutureWarning: remove_na is deprecated and is a private function. Do not use.\n",
      "  stat_data = remove_na(group_data)\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.axes._subplots.AxesSubplot at 0x21410d272b0>"
      ]
     },
     "execution_count": 62,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAdIAAAEKCAYAAACmDdR0AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzs3Xm8lVX9/v/XJc6giIh+HFIcM0VF\nQZwHFM20UNOyshQ1NSvLzMrSykpL036aOYV+DFM/ZpoDjmAEjoCCgIjzQEn6dRbFWbh+f6y14Waz\nz9n7zOfs834+Hudx9r73ute99u2wzn3fa11LtgkhhBBC8yzV0Q0IIYQQurLoSEMIIYQWiI40hBBC\naIHoSEMIIYQWiI40hBBCaIHoSEMIIYQWiI40hBBCaIHoSEMIIYQWiI40hBBCaIGlO7oBoe2tttpq\n7t+/f0c3I4QQupSpU6e+ZrtftXLRkXYD6/RcmTuOOqGjmxFCCO2q33Ffb9H+kv5dS7m4tRtCCCG0\nQJfuSCWtIunb+fVASRMlzZL0iKRDWukYZ+c6z26N+kIIIdSXrn5rdxXg28BFwHvAYbaflrQWMFXS\nGNtvtfAYxwL9bH/YwnpCCCHUoS59RQqcCWwoaTpwtO2nAWy/CLwC9AOQNFvSb/MV6xRJ20gaI+lZ\nSd9qqHJJo4GewGRJh0haQ9KNkmbknx0b2K+/pCckXSbpUUlXSxom6X5JT0saksv1lHS5pIckTZO0\nf2H/eyU9nH92zNt3lzRB0vW5/qslqYE2HJO/65TX573d3PMbQgihiq5+RXoyMMD2wOLG3FEtCzxb\n2PyC7R0knQuMAnYClgdmAZdUqtz2cEnzSvVLuha42/aBknoAvRpp20bAl4BjgIeArwE7A8OBnwEH\nAKcA/7J9pKRVgAcl/ZP0R8Betj+QtDFwDTA417s1sDnwInB//h73VWj7SGAkwMD1NohFZ0MIoY10\n9Y50CZLWBK4EDre9oPDR6Px7JtDL9jvAO5I+kLRKjbeA9wAOA7A9H5jbSNnnbc/MbZoFjLNtSTOB\n/rnM3sBwSSfl98sD65I6yQskDQTmA5sU6n3Q9pxc7/Rc1xIdaQghhPZRVx2ppJWB24BTbU8q+7j0\njHNB4XXpfVuch/JjFI9fOp6Ag2w/WdxR0mnAy8BWpNvvHzRQ73zq7J9hCCF0NV39f8LvACsBSFoW\nuBH4q+3r2uh444DjgPPyrd2etlvyAHIMcLyk4/PV6ta2pwG9gTm2F0g6HOjRkkYv3W/VFs+nCiGE\nUFmXHmxk+3XgfkmPAk8CuwIjJE3PPwMbr6HJvg8Mzbdnp5KeVbbEb4BlgEfyd/hN3n4RcLikSaTb\nuu+28DghhBDaiOwYh1LvtlpvHY89+bsL369x3Mkd2JoQQugaJE21Pbhaua5+axdY+ExxHrAG8AXg\nI9KI3SNaOo9U0qbA3wADB9t+tsouIYQQupEufWu3grtI02G2BJ4CflrLTpK2KNwOLv1Mzh8fANxs\ne+vyTlRS3wr7TZfUt1W/VQghhE6ry16RSjqFNBXlBeBVYKrtsYUik4CDc9kRpA6xBzAA+ANpnuk3\nSKNg9y2fi5r32xc4AZgvaVfbQyUdBpxEukJ9pNJ+ed9RwPvApsB6wBHA4cAOwGTbI3K5vYFfAcux\n6Cp6nqRfkK6uVwAeAI7NA5ImAJOBoaRkp6Ns39uUcxdCCKH1dMkrUkmDgK+Qwgm+CGxbodiRwB2F\n9wNIoQhDgDOA92xvDUwkzw0tZ/t2UljDubkT3ZwUorCH7a1Ig48a04c09/QHwC3AuaQBSlvkbODV\ngFOBYba3AaYAJ+Z9L7C9re0BpM7084V6l7Y9hNTJ/7LSgYvJRm/Mi7FKIYTQVrrqFekuwI2234OF\nUX4L5avVT4CrC5vHF0IY5pI6NkgBDVvWeNw9gOttvwZg+40q5W8phDC8XBbQ0B9YB9iMNPIY0lXy\nxLzvUEk/BlYEViUlMJXafEP+PZVF4Q6LKSYbbbXeOjGiLIQQ2khX7Ugh3VpdQp53+XlgTy8+JLmW\ngIRq1NBxG1AtBGI+cJftry52EGl50hSYwbZfyIOplq9QbwQyhBBCB+uSt3aBe4ADJa0gaSXSs0Qk\n7QP8BBheulptZeOAL5cGE0latYX1TQJ2krRRrm9FSZuwqNN8TVIv8rPeEEIInU+XvJqx/XAOkJ8O\n/BsoDba5gDRo5658q3SS7QZXd2nGcWdJOgO4W9J8YBowogX1vZoHQl0jabm8+VTbT0m6lHTbeTYp\n9L7Zlun3PzF3NIQQ2kgEMnQDgwcP9pQpUzq6GSGE0KXUGsjQVW/thib48JVnePZP+/Psn/bv6KaE\nEELdiY40k3RhDlP4b160e7qkIxopf0J+pnlKLvu2pJn59Slt2M7TSsuuSRohaa22OlYIIYTquuQz\n0rZg+zsAkmaTRvy+VmWXE4CrbJ9BmpfaEUYAj5LWLw0hhNABuvUVqaSekm7LV6CPSvolsBYwXtL4\nXObiHGwwS9Kv8rbvVSg3OwcsIOnEXN+jkk7I2/pLelzSpbmusZJWaKRtR0t6KLftH5JWLPv8YGAw\ncHW+Cm6wrhBCCG2nW3ekwD7Ai7a3yglC55Gu7obaHprLnJIfNm8J7CZpS9vnVygHLExdOgLYDtge\nOFrS1vnjjYELbW8OvAUc1EjbbsjJRlsBjwNHFT+0fT0pCelQ2wNtv1/WjkKy0UdNOyshhBBq1t07\n0pnAMElnSdrF9twKZb4s6WHSVJfNSUlEjdmZlLr0ru15pBSiXfJnz9uenl83mEqUDZB0b05FOpQm\nrn1qe6TtwbYHr9pr2absGkIIoQm69TPSPF9zELAv8DtJxdB7JK1PCqjf1vabOYh++SVrWowa+ayY\nbjSflKHbkFHAAbZn5Lmmu1c5bgghhA7Qra9I84jX92xfBZwDbAO8A6yUi6wMvAvMlbQG8LnC7sVy\nRfcAB+QRvT2BA1kUGNEUKwEvSVqGdEVaSUNtCCGE0E669RUpsAVwtqQFwMfAcaRlzu6Q9FJe8WUa\nKTD+OeD+wr4ji+VKG3Pq0ijgwbzpMtvTJPVvYtt+Tlou7d+kW9CVOsxRwCWS3gd2KH9OWrLc6hux\n4fE3N/HwIYQQahHJRt1AJBuFEELTRbJRWGjeq8/wwMjPVy8YQgihyaIjLVNKLGpqOUm3S1qlGccr\nJSoVf46QNCrPFW1s30g2CiGEDtbdn5FWcgJwFVBtGbbFytnetzkHKyUqlZO0Ww27jyCSjUIIoUN1\n6yvSzpxsVNbOX+SUo0cljVQSyUYhhNAJdOuOlM6dbFR0QU45GkCae/r5piQbvRXJRiGE0Ga6e0fa\nmZONioZKmpxTjvaghpSjYrLRKpFsFEIIbaZbPyPt5MlGpTYsD1wEDLb9gqTTamhDCCGEdtKtr0g7\nebJRSanTfE1SL6A4kjeSjUIIoYN16ytSOneyUam+tyRdSroNPRt4qPDxKGpINurVbyN2PObW5hw+\nhBBCFZFs1A1EslEIITRdJBuFhea+9nRHNyGEEOpWl+xIc6LPBfn1upLGS5om6RFJ++bt/SW9X0gL\nuqSw/yBJMyU9I+l8SQ0OEJK0ad5/mqQN2+C7VEw2au3jhBBCaBv18Iz0VODvti+WtBlwO4umlTxr\ne2CFfS4GjgEm5fL7AHc0UP8BwM22f9mqrc4aSjYKIYTQNXTIFamkUyQ9Kemfkq6RdFID5SZIOk/S\nAznVZ0iFYiaNrgXoTZW4PElrAivbnuj0gPivpM6yUtl9SVGA3ywkGB2Wr3xnSLqykeOMyqlI4yU9\nJ2k3SZfndKNRhXJ7S5oo6WFJ1+WRuRXTjArn5CxJD0p6StIuDTQhhBBCO2j3jjTP2/wKsDXwRWDb\nKrv0tL0j8G3g8gqfnwZ8XdIc0tXl8YXP1s+3ZO8udDhrA3MKZebkbUuwfTtwCXBuHsG7OXAKsIft\nrYDvV2l7H1KAwg+AW4BzSWEKW0gaqBQpeCowzPY2pKSiE/O+S6QZFepd2vYQUidf8Uq5mGw0N5KN\nQgihzXTErd1dSMk/7wFIGl2l/DUAtu+RtLKWXGHlq8Ao23+QtANwpaQBwEvAurZfz533TbkjrPQ8\ntNahy3sA19t+LbfpjSrlb7HtnEj0su2ZAJJmkW4/r0NKSro/X3AuC0zM+w6V9GNgRWBV0hScW/Jn\nN+TfDaYj2R5JmqLDxv17x9DsEEJoIx31jLQp/2MvL1v+/ijSM05sT8xJQKvZfoWcJGR7qqRngU1I\nV6DrFPZfh9pXT1ET215KMlrA4qlGC0jnfj5wl+2vLnaQ6mlGpbrmUx/PuUMIocvqiGek9wAHSlpB\n0krAF6qUPwRA0s7A3Ap5uP8B9sxlPkPqcF6V1E9Sj7x9A1Jg/HO2XwLekbR9fu54GHBzjW0fR8re\n7ZvrXbXG/RoyCdhJ0ka5vhUlbULjaUYhhBA6kXa/msnJP9cC04F/Uz0+701JD5AGFB1Z4fMfApdK\n+gHpanFEvp26K/BrSZ+Qrty+VbgVexwpFWgF0mjdhkbslrd9lqQzgLslzScF2Y+oZd8G6ntV0gjg\nGknL5c2n5gzghtKMmqz3ahu3ZPcQQgiN6PBko3zbcp7tcyp8NgE4yXbE8rRAJBuFEELTKZKNQsnr\nrz/V0U0IIYS61RHTXyZIWtjD2z6NNE1liXQf27s3dDUqaV7h9fzCfqML29dXWsfzaUnXSlo2b18u\nv38mf95fjSQMKc11fSTfPi5vxykV9jul9c5YCCGEzqxTjPhshXSf9xtIMDqLNAf0b0oRgUeRUo2O\nAt60vZGkrwBn2T6kUsWS/gfY0fZ6DbT9DOCMFrY/hBBCF9WmV6SSfi7pCUl3afEEoy/VksyjlKl7\ns6Q7lZKQao7pyyNy9wCuz5uuYFGC0f75PfnzPUvJQRWMBVbPV5q7SNpIKZFphlIaUcX8XUm75yCI\nv+fveaakQ/P3nlnaL48u/odSitFDknbK24coJTpNy78/XTgnN+Rz8rSk39d6TkIIIbS+Nrsizbdv\nDyIlGC0NPEwKEICczKMUwfdLYFgjVQ0BBgDvAQ9Juq3C7d7lJU0BPgHOtH0T0Bd4y/YnuUwxwWht\n4AUA259ImpvLv1bh+MOBW0tXvJIm52PcqDTfs7E/RrYCPgO8QVrP9LL8vb9PSmA6Afgj6ar5Pknr\nAmPyPk8Au+b2DQN+SzqfAANJ5/VD4ElJf7L9QvHAko4h5QnTt29xCmoIIYTW1Ja3dncmhb2/DyDp\nlsJnVZN5Cu6y/Xqu44Zcb3lHuq7tF5Xmi/5LKUno7Qp1lYYoNyvdSGne69q2bwSw/UGVXR7K81ZR\nCoQYm7fPBEqLgQ8DNitcEK+cj9MbuELSxrltyxTqHVeaTyvpMWA98h8GC79MIdlo/fVXjmSjEEJo\nI215a7fBpcloWjJPtWQjbL+Yfz8HTCBdrb0GrCKpVH8xwWgO8CmA/Hlv0lVjNY19p0rK04yKSUel\ndi0F7GB7YP5Z2/Y7wG+A8Tlr9wtUTjaCSDcKIYQO1ZYd6X3AFyQtr5TOs18z69lL0qqSViA947y/\n+KGkPqUwA6UQ+J2Ax/LKLuNZlAp0OIsSjEbn9+TP/+UaJtTafhuYI+mAfLzlJK3YzO9VMhb4buH7\nlAZN9Qb+m1+PaOExQgghtJE260htP0TqsGaQbuVOAcrj/WpxH3AlKQnpHxWej34GmCJpBqnjPNP2\nY/mznwAnSnqG9Az0f/P2/wX65u0nAic3oT3fAL4n6RHgAeB/mvGdir4HDM7Tax4DvpW3/x74naT7\ngR4tOUDfvpu0sIkhhBAa0qbJRpJ62Z6Xr9ruAY6x/XAT9h9BCm7/brWyoWGRbBRCCE2nGpON2vrZ\n2khJm5Ge713RlE40tJ6X33i6o5sQQgh1q007Uttfq6WcpM+SwhOKnrd9IClcvmZK65V+zfZFTdyv\nYhuAZ4B9gdtt/6jCfluQbj0XfWh7u6YcP4QQQtfU4aH1rU1Sf9K8zwGtVN/bQD/bH1Yt3El9aoPe\nfuG55jyeDiGE7qvWW7v1GFp/JrBhTiI6W9KPc5LQDElnAjQhnWg00BOYLOkQSWtIujHvN0PSjg3s\n1z8nOl0m6VFJV0saJun+nEY0JJfrKenynGg0TdL+hf3vzW17uHScnJY0QdL1uf6rpQYTmUIIIbSD\nepx/eDIwwPZASZ8Dfg5sZ/s9LVqI+2pqSCeyPVzSvEKq0bXA3bYPVFo0vFcj7dgI+BIpXegh4Guk\nMInhwM9IU3lOIU29OTLfkn5Q0j+BV4C9bH+QAxmuAUp/FW0NbE6aE3s/abrPfeUHLyYb9Vktko1C\nCKGt1GNHWjQM+Ivt9wBsv9GMdKKiPYDD8n7zaXw6z/O2ZwJImkVKI3JOXeqfy+wNDNeiDOLlgXVJ\nneQFeU7pfKA4f+VB23NyvdNzXUt0pMVko09t0Lu+7t+HEEInUu8dqVgyCam9boXWkmok4CDbTxZ3\nVFrs/GVSVu9SQLGzj1SjEELoROrxGek7wEr59VjgyFL6kKRVW5hONA44Lu/XQ9LKLWzrGOD40nNO\nSVvn7b2Bl2wvIAVAtCiQIYQQQtupu440B9zfL+lRYE9SutKUfBu0dAu1uelE3weG5tuzU0nPKlvi\nN6Qw+kdye3+Tt18EHC5pEum27rstOcgaq27cokaGEEJoWN1NfwlLimSjEEJouu48/SWUeeHNSDYK\nIYS2UtcdqaQvSZolaYHSQuPFz34q6RlJT0o6Ns87Lf5MLpT9nqTHJV1dVkffCvtNl9S3vb5jCCGE\njlXvIz4fBb4I/Lm4Mef/foX0jHMt4J/AJnlKSyXfBj5n+/nixvw8dmDlXUIIIXQHdXFFmpOAHi28\nP0nSabYfL59aku0P/M32h7lzfAYY0kDdlwAbAKMl/UBSL0l/yWlJj0g6qJF2zZN0lqSpOUlpSE4m\nek7S8FymR05geijXd2ze3kvSuJxsNLMs9ehxSZfmq+2xSmu1lh/7GElTJE15/+2Paj+ZIYQQmqQu\nOtJmWBt4ofB+Tt62BNvfIgUkDLV9Likpaa7tLWxvCfyrkeP0BCbYHkSalnM6sBdwIPDrXOaoXN+2\nwLbA0ZLWJ80dPdD2NsBQ4A+FOMCNgQttbw68BSzRmdseaXuw7cErrLxsldMRQgihuer91m5DKoUy\n1Dp8eRjptnDayX6zkbIfAXfm1zNJq8J8XCHdaEtJB+f3vUkd5Rzgt5J2JYU4rA2skcs8b3t6fj21\nUFcIIYR2Vi8d6ScsfnVdLVx2DvCpwvt1SFedtaiUltSQj71oftHCdCPbCyQV042Otz1msYOkRc37\nAYNy5zubRd+rPN1oiVu7IYQQ2ke93Np9GVg9j6JdDvh8lfKjga/kVKP1SVeAD9Z4rLHAd0tvJPVp\nToMLxgDHSVom17eJpJ6kK9NXcic6FFivuQf4VJ8IZAghhLZSFx2p7Y9JzxwnA7cCTwBIOlDSHGAH\n4DZJY3L5WcDfgcdIt16/08iI3XKnA32UlkebQXp+2RKX5XY8nAdM/Zl0p+BqYLCkKcChpe8UQgih\nc4lko26g90b9PPeZVzu6GSGE0KVEslEIIYTQDrpkRypphKQL8ut1JY2XNC3Pw9w3b+8v6f1C2tAl\nhf0H5bmZz0g6X0lDKUXb59/TJG3YQHsmV9hvi/Y5GyGEEDpSPYzaPRX4u+2Lc2LR7SyaDvKs7UrJ\nQxcDxwCTcvl9bN9BhZQiSScDN9v+ZUMNsL1dy75CCCGErqrdr0gbSiFqoOwESedJeiAP7qmUPmSg\ntC5ob6pMY5G0JrCy7Yl5aspfgQMaKLsvcALwTUnj87bD8pXvDElXNnKcUZIuzlfLz0naTdLlOZVo\nVKHc3pIm5gSj6yT1ytt/kdOOHpU0shTGkM/JWZIelPSUpF0aOP7CZKOP3v6gUpEQQgitoCvc2u1p\ne0dS3u3lFT4/Dfh6Hp17O3B84bP18y3ZuwsdztqkeaQljaUa3Q5cApxre6ikzYFTgD1sb0Van7Qx\nfYA9gB8AtwDnkvJ9t5A0UNJqpCvqYTnBaApwYt73Atvb2h5AmidanNKztO0hpE6+4pVyMdlo2ZWr\nTasNIYTQXF3h1u41ALbvkbSypFXKPv8qMMr2HyTtAFwpaQDwErCu7dclDQJuyh1hS1KN9gCut/1a\nbtMbVcrfYts5yehl2zMBJM0i3X5eB9iMtBA5wLLAxLzvUEk/BlYEVgVmkTpjgBvy70g1CiGEDtYR\nHWlTU4jKO7ny90cB+wDYnihpeWA126+wKEloqqRngU1IV6DrFPZvq1QjWJRAtIDF04gWkM79fOAu\n219d7CDpO1wEDLb9Qr71XTxPpbrm0zX+GAohhLrVEbd2m5pCdAiApJ1J4e5zyz7/D7BnLvMZUofz\nqqR+knrk7RuQ0oues/0S8E4ejSvgMODmGts+Dviy8nqjklatcb+GTAJ2krRRrm9FSZuwqNN8LT8z\nPbihCmqx8SrNDkUKIYRQRbtfzeTIu1IK0fNUT+x5U9IDpAFFR1b4/IfApZJ+QLpaHJFvp+4K/FrS\nJ6Qrt28VbsUeB4wiPXu8I//U0vZZks4A7pY0H5gGjKhl3wbqe1UpU/ea/EcFwKm2n5J0KSnofjbw\nUHOPEUIIoW116mQjSROAk2xP6ei2dGW9N1rHc5+ZU71gCCGEhRTJRiGEEELb6xQDVSRdCOxUtvmP\ntndvpfpHAbfavr7CZ2cD+5LWDi0f0ftHUgj+raQRtd+zfW/Z/qcAXyrb7zrbZ7RG20MIIXRuNXWk\neQDMxcAatgdI2hIYbvv01miE7e+0Rj3NdCzQz/aHlT6U9BXgCduHV/o8d5jRaYYQQjdV663dS4Gf\nAh8D2H4E+EpbNapWknpKui2nDD0q6RBJZ0p6LKcPnVMovmtOSHpO0sF5/9FAT2CypEMq1D8Q+D2w\nb87PXUHSPjmFaIakcY207TRJV0gaK2m2pC9K+r1Sxu+dWrT+6KAcGDFV0hil5CUkHZ2TjWZI+oek\nFfP2UUr5wIt9lwrHLyQbvdvMMxxCCKGaWjvSFW2XL3z9SWs3phn2AV60vVVOAJoEHAhsbntL0tqh\nJWsCO5Om25wJYHs48L7tgbavLa/c9nTgF8C1ObO3F+mPioNyslH5Ld1yGwL7AfsDVwHjbW8BvA/s\nlzvTPwEH2x5ESm4qXd3ekJONtgIeJ82XbfC7VGh7IdmoZ5VmhhBCaK5an5G+prTyiQHyVdBLbdaq\n2s0EzpF0Fuk55kTgA+AySbflbSU32V4APCZpjWYeb3vgHtvPQ03JRnfk6T4zgR6kRcRL7e4PfBoY\nANyVk416sOi8DpB0OrAKqQMf08rfJYQQQiuotSP9DjAS2FTSf0nzPw9ts1bVKM+3HEQaLPQ7YCww\nhBTQ8BXgu6RYP1g8WahSTGAtmpVsZHuBpI+9aK5RKdlIwCzbO1TYdxRwgO0Zea7p7uX1FtoUQgih\ng1S9tStpKVJU3TCgH7Cp7Z1t/7vNW1eFpLWA92xfBZwD7Ar0zmHzJ1BhWbQWmgjsJmn9fPyWJhs9\nCfTLGcFIWibnAQOsBLyUb/+26I+WjVf5n5a1MoQQQoOqXpHmq6nvktb87GyjVrYAzpa0gDQQ6kTg\n1pxVK9KqK60mJxEdA9yQ/8B4BdirBfV9lG+Tny+pN+mfx3mkgPqfk9Kf/k26FbxSS9sfQgih9dWU\nbCTp56QBMtcCCzvTGp4Rhk6g94bree6zHX4DIYQQupTWTjY6kvSc9B7S0l1TSWtn1gVJJ5SmlzSl\nnKTbteSybi1px6iGprOEEELonGoabGR7/bZuSAc7gfSscr+y7eUJRSeQprG8B2B7X0lHSCpf4Pv+\nDg6ZCCGE0E5qTTY6rNJ2239t3ea0PUk9gb+T1iHtAVwHrAUMB16zPVTSxcC2wNckLWv7l5K+l8uN\nl1QqN5s0EOsvkk5k0eo0T+dj9SetLHMfsCPwX2B/2+/X0M49SQOoliat/nKc7Q/zMa8AvgAsA3zJ\ndrUVdEIIIbSRWm/tblv42QU4jdTxdEXlIQ7nkRb2Hmp7aC5zSr4vviVplO6Wts+vUA5I6UTAEcB2\npLmmR0vaOn+8MXCh7c2Bt4CDqjUwD5YaBRySAxyWJi39VvKa7W1IsY0nNVBHIdloXrVDhhBCaKaa\nOlLbxxd+jga2JoW4d0UzgWGSzpK0S4WFwiEt3v0wab3RzYHNqtS5M3Cj7XdtzwNuIP3BAfB8TkiC\n9Gy5fw1t/HTe76n8/grS1J6SG6rVt3iyUa8aDhlCCKE5mrv6y3ukK60upzzEQdLY4ud5juhJwLa2\n38wrxyxfpdrGQhGK4QnzSYuJV1MtZKFU53w6yQo+IYTQXdV0RSrpFkmj88+tpCCB0W3btLZRIcRh\nG+AdFs3TXJk0xWdujt/7XGH3Yrmie4ADJK2Yn8EeCNxboVytngD6S9oov/8GcHcL6gshhNBGar2a\nKa6i8gnwb9tz2qA97aE8xOE4YAfgDkkv5UFE00ihCM8B9xf2HVksV9po++F85VoK9r/M9rQ82KjJ\nbH8g6QjgOkmlwUaXNKcugI379GvuriGEEKqoNZDhLNs/qbYtdE6DBw/2lCl1M+03hBDaRWsHMlSK\nwftchW2hE3rmzQigCiGEttJoRyrpuLwE2KeVFsou/TwPPNI+TazYrhGSLii8flVp4e3pkr5ZKHe4\npKfzz+GF7YOUFth+Ji+S3eDgHkmb5nqn5aXkWqP9FxbaW/o5ojXqDiGE0L6qPSP9P1KgwO+Akwvb\n3+lkObvX2v5ucUNemeWXwGDS0mdTJY22/SZp/uUxpIXAbyfNLb2jgboPAG62/cvWamykHoUQQv1o\n9IrU9lzbs21/NS+b9j6pU+olad3mHFBSf0mPFt6fJOm0BspOkHSepAckPSppSBMO9VngLttv5M7z\nLmAfSWsCK9uemNcH/Sups6x0/H1JsYDflDQ+bzssX5XPkHRlI99zlKSLJY2X9Jyk3SRdLunxPDCp\nVG5vSRMlPSzpOkm98vZfSHoof++RpavmfE7OkvSgpKck7dJAE0IIIbSDWqe/fEHS06QFve8GZtPw\nFVxr62l7R+DbwOUNlDkod27XS/pU3rY28EKhzJy8be38unz7EvK6ppcA5+bRvJsDpwB72N4KKM/Y\nLdeHtLD4D4BbgHNJAQ9bSBpqYju4AAAgAElEQVQoaTXgVGBYTiqaQloKDuAC29vm9KUVgM8X6l3a\n9hBSJ1/xSnnxZKO3qzQzhBBCc9U62Oh0UvTdUznAfk8WnxbSlq4BsH0PsLKWXG3lFqC/7S2Bf5JS\ngKByqIEb2V6LPYDrbb+W21Tt9vYt+ap3JvCy7Zm2F5Cm1vQnndPNgPslTQcOB9bL+w6VNDk/o96D\n1AGXNDHZaOUav14IIYSmqrUj/dj268BSkpayPR4Y2MxjflJ23GqpQeWd3GLvbb9uu5T0cykwKL+e\nA3yqUHQdUlbunPy6fHstVKE9jSm1awGLJxwtID2fFun288D8s5nto3LW7kXAwTlr91IWP0+RbBRC\nCJ1ErR3pW/nZ3b3A1ZL+SOoQm+NlYHVJfSUtx+K3LCs5BEDSzsDc8mzc/MyzZDjweH49BthbUh9J\nfYC9gTG2XwLekbR9fu54GHBzjW0fR8rh7ZuPvWqN+zVkErBTKcEoJyNtwqJO87V83mON0hBC6KRq\nvZrZnzTQ6ATgUKA38OvmHND2x5J+DUwmPXOttgTYm5IeIEX3HVnh8+9JGk7q2N8ARuTjvCHpN6RU\nIIBfF27FHkdaXWUF0rPemp732p4l6QzgbknzSaH2I2rZt4H6XpU0Argm/1EBcGrOA76UdEt4duE7\nNMtGfVra34cQQmhITclGAJLWAza2/U9JKwI9bL/Tpo2TJgAn2Y5YnhaIZKMQQmi6Vk02knQ0cD3w\n57xpbeCm5jcvtKdn3nyro5sQQgh1q9Zbu98BhpBux2L7aUmrt1YjJF0I7FS2+Y+2d2+tY+TjnACM\ntP1erW0AniFNgfkY2MH2+2X7nQJ8qWy/62yf0WoNDyGE0GnVGlo/2fZ2kqbZ3lppRZKH85STLkPS\nbGBwafpK2Wc9bM+vsP0SYLLtv7RDE9vEKhtu5LeefaajmxFCCF1Kq97aJQ2u+RmwgqS9gOtI8zc7\nLUk9Jd2WE4gelfRLYC1gfCGlaJ6kX0uaTFpKrbyObwJfBn4h6eq87cdKOb0zJJ3ZyPEnSDpX0j05\nzWhbSTco5f6eXij39ZxSNF3SnyX1yNsvzoEKsyT9qlB+tqRf5SSkmZI2baVTFkIIoRlqvbV7MnAU\naRTpsaR82svaqlGtZB/gRdv7AUjqDRwBDC1ckfYEHrX9i0oV2L4sT7u51fb1kj5HihPczvZ7NUx/\n+cj2rpK+T5piM4g0svhZSecCq5Om9+yURzNfRBoV/VfglDzyuAcwTtKWtksLBbxmextJ3wZOAr5Z\nfmBJx5DyhFlhtViPNIQQ2kqjHamkdW3/J6fxXJp/uoqZwDmSziJ1hPdqyUVe5gP/aEKdw4C/lJ6x\n1pBsNLrQlll5DiuSniOFRexM6lwfym1bAXgl7/Pl3BkuDaxJSkAqdaTFZKMvVjqw7ZGkhchZZcON\nmhIiEUIIoQmqXZHeBGwDIOkftg9q+ya1jjwXcxCwL/A7SWMrFPug0nPRRrRFstEVtn+62EGk9UlX\nmtvaflMp5D6SjUIIoROq9oy0eAm3QVs2pLVJWgt4z/ZVwDmkPwjeAVZqQbVjgSPzPNrWSDYaBxxc\nGgEtadU8X3dl4F1grqQ1iEXUQwih06p2NeMGXncFWwBnS1pAmrpyHGlA0R2SXrI9tKkV2r5T0kBg\niqSPSM+Kf9bcBtp+TNKpwFhJS+V2fsf2JEnTSOH2z9HCBQI26lOe8x9CCKG1NDr9JcfgvUu6Ml0B\nKM2/FGDbsaxIFxDJRiGE0HS1Tn9p9IrUdo/Wa1LoKM++2aZJjiGE0K3VOo+0U8lzKVcrvJ6Z52FO\nKZRZVdJded7mXXkFGJScL+kZpcXAS4Opbsx1FH8+K+nsPJfz7AbacmGF/Y5oj/MQQgih49XLiM/i\n3NCSk4Fxts+UdHJ+/xPSwJ2N8892wMWkeaEHVqpY0nVAv8Kap4ux/Z1W+g4hhBC6oHa/IpV0Vg4S\nKL0/TdIPK5TbPacC3SjpMUmX5AE5tdofuCK/voIUpFDa/lcnk4BVtPiapsU2jCaFNkyWdIikNXJ7\nZuSfHRvYr7+kJyRdppSqdLWkYZLuz1fIQ3K5npIul/SQpGmS9i/sf29OL3q4dJx8TiZIuj7Xf7Uq\nTI4NIYTQfjri1u7fyIt1Z18mRQ5WMgT4IWkE7oZUDh8wadTr1BxgULJGKQAh/y6F7K8NvFAoNydv\nW7Jiezjwvu2Btq8Fzgfutr0VaTrNrAa/JWxECr3fEtgU+BopgOEkFo30PQX4l+1tgaGkUcY9SaEM\ne9nehnSuzi/UuzVpXdjNSFOSyoP2gZRspBQxOOXDt+dWKhJCCKEVtPutXdvTJK2e53n2A960/Z8G\nij9o+zkASdeQOqLry8rsZPvFPBfzLklP2L6nkSZUuoKrdWrPHsBh+XvMBxrroZ63PTO3fRbpNrMl\nzQT65zJ7A8MlnZTfLw+sC7wIXJCn2swHNinU+6DtObne6bmu+5b4QoVkoz4bfrqrTV0KIYQuo6Oe\nkV4PHAz8D+kKtSHlHcASHYLtF/PvVyTdSLqKvQd4WdKatl/Kt25L0XtzSPF8JeuQOq7WVp5kVEw5\nKp13AQfZfrK4o6TTgJeBrUh3DT5ooN5INgohhA7WUaN2/wZ8hdSZll9hFg2RtH5+NnoIZVde+Rnj\nSqXXpCu8R/PHo4HD8+vDSaHxpe2H5dG72wNzS7eAazCOFOyApB6SWjqPdgxwfOk5p6St8/bewEs5\n4/gbQExDCiGETqpDOlLbs0hRff+t0olNBM4kdY7PAzeWfb4GcJ+kGcCDwG2278yfnQnsJelpYK/8\nHlIa0XOkBbsvBb5N7b4PDM23Z6cCmzdh30p+AywDPCLp0fwe4CLgcEmTSLd1323JQTbs05JUxBBC\nCI2paWHvjiBpd+Ak25/v6LZ0dZFsFEIITVdrslGXDGQITfPcm+93dBNCCKFudfhAFUlbAFeWbf7Q\n9nbAhAb2mQ0Mtv2apB+QFrY2ad3PI2x/oLT02G4sGlk7wvb0/Dzyj6Tl1d4DRpDC4iu2ISca7Qvc\nbvtHZe3oS3puWm5P269X++4hhBC6vg7vSPMUkYHN2VfS2sD3gM1svy/p76RBTKNykR/ZLh/MtESy\nUe60G2rDsTSQbJQ7y2a1PYQQQn2oh2SjpYEVJC0NrEj1qSyRbBRCCKHVdOlkI9v/JS3a/R/gJdJU\nlrGFImcoBdOfK2m5vK0bJhu92UgzQwghtES7d6S2pwGrS1pL0lbUkGyUU4RKyUYLKa3osj+wPrAW\n0FPS1/PHPyV1YNsCq5IC66HlyUYX5+8x33bVZKM8F3RhshHpOW7/XGZv4OScUDSBRclGywCX5mk2\n15E6zZIHbc/J9U4v1LX4F7JH2h5se/ByK/ep8euFEEJoqq6ebDSM1GG9CiDpBmBH4KrC/NQPJf2F\ndCUIkWwUQgihFXXpZCPSLd3tJa2YnxXuCTwOUHrumbcfwOKJR5FsFEIIoVV06WQj25NJHfHDpFum\nS5GD2oGr863RmcBqwOl5e7dLNtqgzwotamQIIYSGRbJRNxDJRiGE0HSRbBQWeumtjzu6CSGEULc6\nfKBKc5KNmlD3ncD2wH3FK1tJ65Oe065Kui18JvCXSm2Q1A+4FVgW+J7tewv1RLJRCCF0cx3ekbYk\n2agGZ5NCGo4t234WcK7tv0m6BNjedkNt2BN4wvbh5R9EslEIIYS6uLXbUFqS7XHAO2VlRZoPWhot\nfAVpVG+legcCvwf2lTRd0gqS9slpQzMkVboaLbbhCkljJc2W9EVJv5c0U9KdkpbJ5QZJulvSVElj\nCqONj86JRzMk/UPSinn7KEnnS3pA0nOSDm72iQshhNBiddGR0rS0pL7AW7Y/ye8bSzaaDvwCuDZf\nsfYijfQ9KKcbfalKuzYE9iOFRlwFjLe9BfA+sF/uTP8EHGx7EHA5cEbe9wbb2+bjPA4cVah3TVI4\nxedZtM7qYorJRu++HXeZQwihrXT4rd3WYHuapNUlrQX0o/G0pJYkG20P3GP7+XzcN6qUv8P2x3m6\nTA+gtOh4Kd3o08AA4K48lbQHKeoQYICk04FVSB34mEK9N+U5po9JWqPSgW2PJE8FWnujrTrn0OwQ\nQqgDddGRZrWmJb1GCqpfOl+VNiXZSNTe6UJOIbK9QNLHXjTXqJRuJGCW7R0q7DsKOMD2DEkjgN3L\n6y20KYQQQgepl1u7UGNaUu7MxudyAIcDN9d4jInAbnnUL5JWbXZrkyeBfpJ2yPUtI6kU8rAS8FK+\n/XtoC48TQgihjdRNR1opLUnSvaRnpXtKmiPps7n4T4ATJT1Demb6vzUe41XgGOAGSTOAa1vY5o9I\nHfpZub7ppKxggJ8Dk4G7gCdacpw1V1mmJbuHEEJoRKdNNgqtZ8MNB/rZZ6d3dDNCCKFLiWSjEEII\noR3U02CjFpF0CktOZ7mONA/1OOBh20s8q5R0BCnMvuh+299pk4aGEELoVOLWbhWSngA+V5ry0hXF\nrd0QQmi6ur21K+kwSY/kxJ8rJX1B0mRJ0yT9szSvMicLXSnpX5KelnR0oY4f5dSgRyT9qpFjXQJs\nAIyW9ANJvST9JacTPSLpoEb2nZcTl6bmdg2RNCGnEQ3PZXpIOrvQlmPz9l6SxuUEpZmS9s/b+0t6\nXNKlkmbl1KSKa6QVAxnejkCGEEJoM13qijRPDbkB2Mn2a3n6iUlJRZb0TeAztn8o6TTgQFKIQk9g\nGrAdKQDhYFL+rkgLff/e9j0NHHM2MDgf7yxgOdsn5M/62H6zgf0M7Gv7Dkk35jbsB2wGXGF7oKRj\ngNVtny5pOeB+0u3lF4AVbb8taTVgErAxsB5pHdXBtqdL+jsw2vZVjZ23uCINIYSmq/WKtKs9I90D\nuN72a5CShZRWj7k2Z9QuS1oAvORm2+8D70saDwwhRevtTepYIaUGbQxU7EjLDCPNVSUfv2Inmn3E\n4klGHxZSjvrn7XsDWxbycnvntswBfitpV1J4w9pAKcHo+RxdCGlx8VJdIYQQOkBX60grJQv9Cfj/\nbI9WWgz8tMJn5WWd6/id7T+30vEbUp5kVEw5Kp13AcfbLsb/kZOM+gGDcuc7G1g+f1xMNZoPVLy1\nG0IIoX10tWek44AvK60DWkoW6g38N39evtTZ/pKWz+V3Bx4iZdYeKalXrmNtSavXePyxwHdLbyT1\nae4XycYAx2nRSjCbSOpJ+k6v5E50KOmWbrP16dPV/l4KIYSuo0v9H9b2LElnAHdLmk+6PXsacJ2k\n/5KeJa5f2OVB4DZgXeA3tl8EXpT0GWBiDoqfB3wdeKWGJpwOXCjpUdLV4K9Iz2yb6zLSrdmHlRrz\nKmlJt6uBWyRNIaUdtSjZKIQQQtvpUoONmiIPNppn+5yObktH23T9gX7i+RhsFEIITVG3019CCCGE\nzqRL3dptCtun1Vo2P0MdV+GjPW03OglT0mRgubLN37A9s9bjl9Un0p2CBc3ZP4QQQvuq2460KXJn\nObCZ+27X1H0knQgcmd9eBtwE3EFa3m0H4ABJOwI/I43svc32T/K+R5FWr3kReJo0rea7hBBC6BBx\na7edSRoEHEEKh9geOBroA3wa+KvtrYGPgbNI82YHAttKOkDSWqTl1bYH9gI2beQ4C5ON3nonko1C\nCKGtREfa/nYGbrT9ru15pFG/uwD/tj0pl9kWmGD7VdufkEbx7koKlLjb9hu2PyaF6ldke6TtwbYH\nr7JS3zb9QiGE0J1FR9r+1MD2d2so09D2EEIIHSQ60vZ3D+kZ6Io5fOFA4N6yMpOB3SStJqkH8FXg\nbtK82N0k9cnpSA2G5ocQQmgfMdiondl+WNIoUqcIabDRm2VlXpL0U9LgIwG3274ZQNJvSR3ti8Bj\nwNxqx+zVN/4xhxBCW6nbQIZ6JamX7Xn5ivRG4HLbNza2z+DBgz1lypT2aWAIIdSJCGSoX6dJmg48\nSlrp5qZqO3zw6sdt3qgQQuiu4p5fFXmd0NuA1Uirxlzbke2xfVJHHj+EEMLioiPNJC2dp5qU2xpY\nxnazAhtCCCHUt7rsSPNo2L8D6wA9gN+QOsThwCfAWNsn5UE/b+TPHgZ+WFbP6sBVQL98O/UgYFXg\nj0BP0tqge9p+p0IbRpBWcukBDAD+QFp4/Bt5v33zwuQbAheS1h99Dzja9hOSvgCcmvd5HTjU9ss5\njH9dYIP8+zzb57fwlIUQQmimuuxIgX2AF23vByBpPeDXwKa2LWmVQtlNgGG255dXYvsVSd8ETrL9\neUnLAncBh9h+SNLKwPuNtGMAqZNeHngG+IntrSWdCxwGnAeMBL5l+2lJ2wEXkRKN7gO2z+39JvBj\nFnX0mwJDgZWAJyVdnAMaFpJ0DHAMwFqrrlPTSQshhNB09dqRzgTOkXQWcCswEfgAuEzSbXlbyXWV\nOtEGfBp4yfZDALbfrlJ+fL5afUfSXOCWQvu2zIuL70haT7W0TykAfx3gWklrkq5Kny/Ue5vtD4EP\nJb0CrAHMKR7Y9khSJ82A9baKodkhhNBG6nLUru2ngEGkDut3pPD3IcA/SLdb7ywUf3eJChomoCmd\n0oeF1wsK7xeQ/ohZCnjL9sDCz2dymT8BF9jeAjiWdFVbqd751O8fRCGE0OnVZUeaw93fs30VcA4p\np7a37duBE2jmSi/AE8BakrbNx1kpz+dslnxF+7ykL+X6JGmr/HFv4L/59eHNPUYIIYS2Va9XMlsA\nZ0taQFpJ5UTgVknLk64qf9CcSm1/JOkQ4E+SViA9Hx0GzGtBWw8FLpZ0KrAM8DdgBnAa6Zbvf4FJ\nwPrNPcDy/ZZpQfNCCCE0JpKNuoFINgohhKaLZKOw0Mcvf1i9UAghhGap11u7SJpAmrYyRdI8270q\nlOlHGsG7LHAHsG9ZkfuBCaSpM//P9tAKdXyWtAh30fO2D2zxlwghhNDp1W1HWqM9gSdslwbz/Ky8\ngKQ7gW/bHl+pAttjgDFt18QQQgidWV3c2pXUU9JtkmZIejQPCCovc0b+fJKkNSQNBH4P7Ctpeh48\nVL7PL4CdgUsknS2ph6RzJM2U9Iik4xtp02xJv5U0UdIUSdtIGiPpWUnfKpT7kaSHcn2/Kmy/SdJU\nSbNyuEJp+7zy79L8MxdCCKGl6qIjZVGS0Va2B7D4PFFIcX6TbG9FWlj7aNvTgV8A1+b5m0skFNn+\nNTCFFM/3I1JS0PrA1ra3BK6u0q4XbO9AWrh7FHAwsD3pVjGS9gY2Js1xHQgMkrRr3vdI24OAwcD3\nJPVt6LtUOrCkY3IHPuX1eW9UaWYIIYTmqpeOdCYwTNJZknaxXb7Y9UcsSjOaCvRv5nGGAZeUwu1t\nV+uhRhfaN9n2O7ZfBT7IMYV7559ppKzfTUkdK6TOcwZp6sunCttr+i62R9oebHtw316rNu1bhhBC\nqFldPCO1/ZSkQaTBQr+TNLasyMdeNM+nJUlAzU02KqYald4vnev7ne0/L3YQaXdSp72D7ffywKlS\nslFrfZcQQgitoC6uSCskGW3TRocaC3yrlGYkqaWXemOAI3PmLpLWzivO9AbezJ3opqTbwSGEEDqh\nuuhISUlGD+alzk4BTm+j41wG/Ad4JN92/VpLKrM9Fvg/YKKkmcD1pBVd7gSWlvQIaQm4SS05zjJr\nLFe9UAghhGaJZKNuIJKNQgih6SLZKCz08ctNWeAmhBBCU8RAlUzSZBatBVryDeD/sSj96Hu27y3b\n70aWDJT/SQ5qCCGEUOeiI81sb1dpu6SvsHj6Ufl+EQUYQgjdWN3f2q2UepRTh86S9GD+2aiBfZdI\nP5K0j6SHc33jGjnuaZKukDQ2H++Lkn6fU5HulLRMLjdI0t05xWiMpDXz9qNz4tEMSf+QtGLePkrS\n+ZIekPScpINb/6yFEEKoVd13pDScevS27SHABcB5lXYsTz8CegGXAgflZKEvVTn2hsB+wP7AVcB4\n21uQ1jHdL3emfwIOzilGlwNn5H1vsL1tPs7jwFGFetckRRd+Hjiz0oGLyUZvvPtWlWaGEEJoru7Q\nkTaUenRN4fcONda1PXCP7eehpmSjO2x/nNvQg0Wd+ExSItGngQHAXXnqzqnAOrnMAEn35mkxhwKb\nF+q9yfYC248BFbN2i8lGq/ZcpcavF0IIoanq/hlpI6lHxXk/tc4Balayke0FkoqJRMVko1k5j7fc\nKOAA2zMkjQB2L6+30KYQQggdpO6vSBtJPTqk8HtijdVNBHaTtH6uu6XJRk8C/STtkOtbRlLpynMl\n4KV8+/fQFh4nhBBCG6n7K1JS6tHZkhYAHwPHkRKElstTXpYCvlpLRbZfzUua3SBpKeAVYK/mNsz2\nR3mw0PmSepP+eZwHzAJ+DkwG/k26FbxSc4+zzBo9m7trCCGEKrplspGk2cBg2691dFvaQyQbhRBC\n00WyUVjok1fe7ugmhBBC3eqWHant/uVXo5JOyXNFiz+n5M/OljRL0tnldUk6osJ+F7bXdwkhhNCx\nuuWt3aaS9DbQz/aHVQt3QgPX3cTT//NURzcjhBC6lG59a7eRNKPfSpqYgwq2yUlCz0r6ViN1jQZ6\nApNzPWtIujHXPUPSjg3s11/SE5Iuy224WtIwSfdLelrSkEJbL88pRtMk7V/Y/96covRw6TiSdpc0\nQdL1uf6rJcUUmBBC6CD1Omq3lGa0H0AeEXsW8ILtHSSdS5qnuROwPGmU7CWVKrI9XNK8nGyEpGuB\nu20fKKkHKe2oIRuR0o+OAR4irV+6MzAc+BlwAGn91H/ZPlLSKqR1Vf9JHhFs+wNJG5OCI0p/GW1N\nCmh4Ebg/f4/7igfOo4uPAVinz+rVz1gIIYRmqcsrUhpOMxpd+Hyy7Xdsvwp8kDuxWuwBXAxge36h\n7kqetz3T9gJSZz0uhzKUko0A9gZOzslGE0gd+7rAMsClOdnoOmCzQr0P2p6T651eqGuhYrJR3169\na/xqIYQQmqour0gbSTMqPeNcwOLpQKWkodZWfozi8UvHEym798nijpJOA14GtiL9wfNBA/XOp07/\nOYYQQldQl1ekjaQZtYZxpFAHJPWQtHIL6xsDHF96zilp67y9N/BSvur8BimrN4QQQidTlx0pKc3o\nwXy79BTg9Fas+/vA0HzLdSqLh8k3x29It3EfkfRofg9wEXC4pEnAJsC7zT3A0qu3tK8PIYTQkJj+\n0g1IeoeU6xuS1YBukWpVozgfi4vzsbjufD7Ws92vWqF4ttY9PFnLXKjuQtKUOB+LxPlYXJyPxcX5\nqC460kzSFsCVZZs/tL1dlf36kp6bltvT9uut1b4QQgidU3Skme2ZwMBm7Pd6c/YLIYRQH+p1sFFY\n3MiObkAnE+djcXE+FhfnY3FxPqqIwUYhhBBCC8QVaQghhNAC0ZHWEUn7SHpS0jOSTq7w+XKSrs2f\nT5bUv/1b2X5qOB8nSnpM0iOSxklaryPa2V6qnY9CuYMlWVJdj9Ss5XxI+nL+d2SWpP9r7za2pxr+\ne1lX0vi8uMYjkvbtiHZ2Srbjpw5+SMlHzwIbAMsCM4DNysp8G7gkv/4KcG1Ht7uDz8dQYMX8+rju\nfj5yuZWAe4BJwOCObncH//uxMTAN6JPfr97R7e7g8zESOC6/3gyY3dHt7iw/cUVaP4YAz9h+zvZH\nwN+A/cvK7A9ckV9fD+xZx0uwVT0ftsfbfi+/nQSs085tbE+1/PsBKVnr9yye7VyPajkfRwMX2n4T\nwPYr7dzG9lTL+TBQiknrTVp9KhC3duvJ2sALhfdz8raKZWx/AswF+rZL69pfLeej6CjgjjZtUceq\nej5yzvOnbN/ang3rILX8+7EJsEleQ3iSpH3arXXtr5bzcRrwdUlzgNuB49unaZ1fzCOtH5WuLMuH\nZNdSpl7U/F0lfZ201utubdqijtXo+ZC0FHAuMKK9GtTBavn3Y2nS7d3dSXcr7pU0wPZbbdy2jlDL\n+fgqMMr2HyTtAFyZz8eCtm9e5xZXpPVjDvCpwvt1WPLWy8IykpYm3Z55o11a1/5qOR9IGkZa2GC4\n7Q/LP68j1c7HSsAAYIKk2cD2wOg6HnBU638vN9v+2PbzpLzqjdupfe2tlvNxFPB3ANsTSWsnr9Yu\nrevkoiOtHw8BG0taX9KypMFEo8vKjAYOz68PBv7lPHKgDlU9H/lW5p9JnWg9P/+CKufD9lzbq9nu\nb7s/6ZnxcNtTOqa5ba6W/15uIg1IQ9JqpFu9z7VrK9tPLefjP8CeAJI+Q+pIX23XVnZS0ZHWifzM\n87uk9U0fB/5ue5akX0sanov9L9BX0jPAiUCDUyC6uhrPx9lAL+A6SdMllf+Po27UeD66jRrPxxjg\ndUmPAeOBH7lO87NrPB8/BI6WNAO4BhhRx3+IN0kkG4UQQggtEFekIYQQQgtERxpCCCG0QHSkIYQQ\nQgtERxpCCCG0QHSkIYQQQgtERxpCFyZpfp66U/rp34w6VpH07dZv3cL6hze22kwbHfMASZu15zFD\n9xXTX0LowiTNs92rhXX0B261PaCJ+/WwPb8lx24LObXrMtJ3ur6j2xPqX1yRhlBnJPWQdLakh/K6\nkcfm7b3yuqsPS5opqbS6x5nAhvmK9mxJu0u6tVDfBZJG5NezJf1C0n3AlyRtKOlOSVMl3Stp0wrt\nGSHpgvx6lKSL87qWz0naTdLlkh6XNKqwzzxJf8htHSepX94+MAfIPyLpRkl98vYJkn4r6W7gJ8Bw\n4Oz8nTaUdHQ+HzMk/UPSioX2nC/pgdyegwtt+HE+TzOk/7+9uwu1MYvjOP79aXRGkpea5MqUSEOc\nMlLexlvKjeFibkgKF+ZCKZQblCTlUrlgamZEEUm5mUMuTo2i8XJwSM2FRGlceO2gyH8u1n8327Gd\nc+Y8aM6e3+dmP3s9a+31rOdi/1t7PXv9tSfLeh2v/f9403qzgW2IpI48vhMRyyl7oj6NiOmSWoDz\nks5Qsnssj4hnueXdhdzNaSswOSJaASTN66XPVxExO+ueA9ZHxJ+SZgD7gQW9tB+ZdZYCp4FZwDrg\nD0mtEdEBDAWuRMQmSS9jseEAAAIwSURBVNuBHZSddw4BGyKiXdLOLN+YnzsiIr7L6xpP3YxU0pOI\nOJjHu/Ie7ct2Y4DZwETKtngnJC0BlgEzIuKFpFFZ90A/xmtNzoHUbGB7WQuAdRYDU+pmV8Mpm63f\nB3ZLmgu8paTJGt2PPo9BmeECMylbLNbOtfSh/emICEk3gL8i4kZ+3k3ga6Ajr+9Y1j8MnJQ0nBIs\n27P8V+B49+v6gMkZQEdQtoVsqzt3KjOY3JJUux+LgJ9r+Woj4lGF8VqTcyA1az6izNra3iksP89+\nBUyLiNcqWV6+bND+De8u+3Sv05Wvg4AnDQJ5b2pZdt7WHdfef+g7qS8Pc3T1cO4XYFlEXMv7MK/B\n9cA/6cTUoM/+jteanNdIzZpPG/CjpMEAkiZIGkqZmT7MIDofGJv1n1PSqNXcBb6R1JKzwIWNOomI\nZ8AdST9kP5I09SONYRAlQxHACuD3iHgKPJY0J8tXAe2NGvP+mIYBD/KerOxD/2eANXVrqaM+8Xht\nAHMgNWs+PwG3gCuSOimp4r4AjgDfSrpECSa3ATKjyXlJnZL2RsQ9St7J69nmag99rQTWqmQEuQl8\n30Pdf6MLmCTpMmUNcmeWr6Y8RHQdaK0r7+4osEXSVUnjgG3AReAsOe6eRMRvlPXSS7kGvTlPfarx\n2gDmv7+Y2X+OPsLfesw+F89IzczMKvCM1MzMrALPSM3MzCpwIDUzM6vAgdTMzKwCB1IzM7MKHEjN\nzMwqcCA1MzOr4G/pi1bxaOz3LwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x21410d3e7f0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "sns.barplot(data=fi_df, x='Feature importance', y='Feature')\n",
    "#plt.xscale('log')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As expected, by far the most important feature is `t2m_fc_mean`. Does the ordering of the remaining features make meteorological sense?\n",
    "\n",
    "There some other really cool ways of interpreting a random forest. For more information check out this notebook: https://github.com/fastai/fastai/blob/master/courses/ml1/lesson2-rf_interpretation.ipynb"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Your turn\n",
    "\n",
    "Again, here are some challenges for all of you that are bored ;)\n",
    "\n",
    "1. The different stations have very different characteristics. Try building a model (e.g. linear regression) for each station individually. See if you can improve the validation and test score.\n",
    "2. Try adding some time information (for example the month as a feature. Does this improve the forecast?\n",
    "3. Try to get as good a score as possible. Be creative and use Google!"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.3"
  },
  "toc": {
   "nav_menu": {},
   "number_sections": true,
   "sideBar": false,
   "skip_h1_title": true,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}