Random Forests.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 1. Ensemble of Random Forests with Resampling\n",
    "\n",
    "In this notebook we use an enseble of random forest to predict the customer satisfaction. One of the more poignant characteristics of the data is the imbalance between the target classes, only ~4.5% of the observations are labeled as unsatisfied. To deal with this issue we undersample the majority class. The procedure is as follows (for details see train_forests.py):\n",
    "\n",
    "- Split the training data according to the target variable: happy and unhappy customers.\n",
    "\n",
    "- Randomly pick a specified fraction of the unhappy data.\n",
    "\n",
    "- Merge the previous observations with a equally sized random sample of happy observations.\n",
    "\n",
    "This gives us a balanced data set on which to train our random forests.\n",
    "\n",
    "Remark: The metric that we use for evaluation is the area under the ROC curve."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import pickle\n",
    "import seaborn as sns\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from sklearn.ensemble import RandomForestClassifier\n",
    "from sklearn.model_selection import train_test_split\n",
    "from sklearn.metrics import roc_auc_score"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# We import custom utility functions for data processing and random forest training\n",
    "from process_data import process, create_submission, drop_columns\n",
    "from train_forests import trainForests, mean_ensemble"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#data = pd.read_csv('data/train_saldo.csv')\n",
    "#data = pd.read_csv('data/train_extended_saldo.csv')\n",
    "data = pd.read_csv('data/train.csv')\n",
    "process(data)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "train, test = train_test_split(data, test_size = 0.2, random_state = 42)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "X_train, Y_train = train.ix[:,:-1], train['TARGET']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "X_test, Y_test = test.ix[:,:-1], test['TARGET']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>var3</th>\n",
       "      <th>var15</th>\n",
       "      <th>imp_ent_var16_ult1</th>\n",
       "      <th>imp_op_var39_comer_ult1</th>\n",
       "      <th>imp_op_var39_comer_ult3</th>\n",
       "      <th>imp_op_var40_comer_ult1</th>\n",
       "      <th>imp_op_var40_comer_ult3</th>\n",
       "      <th>imp_op_var40_efect_ult1</th>\n",
       "      <th>imp_op_var40_efect_ult3</th>\n",
       "      <th>imp_op_var40_ult1</th>\n",
       "      <th>...</th>\n",
       "      <th>saldo_medio_var33_ult1</th>\n",
       "      <th>saldo_medio_var33_ult3</th>\n",
       "      <th>saldo_medio_var44_hace2</th>\n",
       "      <th>saldo_medio_var44_hace3</th>\n",
       "      <th>saldo_medio_var44_ult1</th>\n",
       "      <th>saldo_medio_var44_ult3</th>\n",
       "      <th>var38</th>\n",
       "      <th>SumZeros</th>\n",
       "      <th>NumAssets</th>\n",
       "      <th>TARGET</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0</td>\n",
       "      <td>23</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>10.576564</td>\n",
       "      <td>315</td>\n",
       "      <td>4</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>0</td>\n",
       "      <td>34</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>10.805234</td>\n",
       "      <td>289</td>\n",
       "      <td>9</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>0</td>\n",
       "      <td>23</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>11.117417</td>\n",
       "      <td>300</td>\n",
       "      <td>6</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>0</td>\n",
       "      <td>37</td>\n",
       "      <td>0.0</td>\n",
       "      <td>195.0</td>\n",
       "      <td>195.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>11.066763</td>\n",
       "      <td>270</td>\n",
       "      <td>12</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>0</td>\n",
       "      <td>39</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>11.672584</td>\n",
       "      <td>279</td>\n",
       "      <td>11</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>5 rows × 331 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "  var3  var15  imp_ent_var16_ult1  imp_op_var39_comer_ult1  \\\n",
       "0    0     23                 0.0                      0.0   \n",
       "1    0     34                 0.0                      0.0   \n",
       "2    0     23                 0.0                      0.0   \n",
       "3    0     37                 0.0                    195.0   \n",
       "4    0     39                 0.0                      0.0   \n",
       "\n",
       "   imp_op_var39_comer_ult3  imp_op_var40_comer_ult1  imp_op_var40_comer_ult3  \\\n",
       "0                      0.0                      0.0                      0.0   \n",
       "1                      0.0                      0.0                      0.0   \n",
       "2                      0.0                      0.0                      0.0   \n",
       "3                    195.0                      0.0                      0.0   \n",
       "4                      0.0                      0.0                      0.0   \n",
       "\n",
       "   imp_op_var40_efect_ult1  imp_op_var40_efect_ult3  imp_op_var40_ult1  \\\n",
       "0                      0.0                      0.0                0.0   \n",
       "1                      0.0                      0.0                0.0   \n",
       "2                      0.0                      0.0                0.0   \n",
       "3                      0.0                      0.0                0.0   \n",
       "4                      0.0                      0.0                0.0   \n",
       "\n",
       "    ...    saldo_medio_var33_ult1  saldo_medio_var33_ult3  \\\n",
       "0   ...                       0.0                     0.0   \n",
       "1   ...                       0.0                     0.0   \n",
       "2   ...                       0.0                     0.0   \n",
       "3   ...                       0.0                     0.0   \n",
       "4   ...                       0.0                     0.0   \n",
       "\n",
       "   saldo_medio_var44_hace2  saldo_medio_var44_hace3  saldo_medio_var44_ult1  \\\n",
       "0                      0.0                      0.0                     0.0   \n",
       "1                      0.0                      0.0                     0.0   \n",
       "2                      0.0                      0.0                     0.0   \n",
       "3                      0.0                      0.0                     0.0   \n",
       "4                      0.0                      0.0                     0.0   \n",
       "\n",
       "   saldo_medio_var44_ult3      var38  SumZeros  NumAssets  TARGET  \n",
       "0                     0.0  10.576564       315          4       0  \n",
       "1                     0.0  10.805234       289          9       0  \n",
       "2                     0.0  11.117417       300          6       0  \n",
       "3                     0.0  11.066763       270         12       0  \n",
       "4                     0.0  11.672584       279         11       0  \n",
       "\n",
       "[5 rows x 331 columns]"
      ]
     },
     "execution_count": 34,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "data.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1.2 Single Forest Classifier\n",
    "\n",
    "To have a better understanding of the problem, we examine the performance of a single random forest classifier."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "rf = RandomForestClassifier(n_estimators = 150)\n",
    "rf.fit(X_train, Y_train)\n",
    "rf_bal = RandomForestClassifier(n_estimators = 150, class_weight='balanced')\n",
    "rf_bal.fit(X_train, Y_train)\n",
    "# custom undersambpling (see train_forests.py)\n",
    "rf_custom = trainForests(train, n_trees = 150)[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "rf --> Test score = 0.7628327185552672, Train score = 0.9975213935333158\n",
      "rf_bal --> Test score = 0.7553627220686608, Train score = 0.9912921279519079\n",
      "rf_custom --> Test score = 0.8040550522725947, Train score = 0.9400009333039534\n"
     ]
    }
   ],
   "source": [
    "# rf\n",
    "Y_prob = pd.DataFrame(rf.predict_proba(X_test))[1]\n",
    "Y_prob_train = pd.DataFrame(rf.predict_proba(X_train))[1]\n",
    "print(\"rf --> Test score = {}, Train score = {}\".format(roc_auc_score(Y_test, Y_prob), roc_auc_score(Y_train, Y_prob_train)))\n",
    "# rf_bal\n",
    "Y_prob = pd.DataFrame(rf_bal.predict_proba(X_test))[1]\n",
    "Y_prob_train = pd.DataFrame(rf_bal.predict_proba(X_train))[1]\n",
    "print(\"rf_bal --> Test score = {}, Train score = {}\".format(roc_auc_score(Y_test, Y_prob), roc_auc_score(Y_train, Y_prob_train)))\n",
    "# rf_custom\n",
    "Y_prob = pd.DataFrame(rf_custom.predict_proba(X_test))[1]\n",
    "Y_prob_train = pd.DataFrame(rf_custom.predict_proba(X_train))[1]\n",
    "print(\"rf_custom --> Test score = {}, Train score = {}\".format(roc_auc_score(Y_test, Y_prob), roc_auc_score(Y_train, Y_prob_train)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "From this we see that the built-in parameter to deal with imbalanced data sets, 'class_weight', is not really useful. However, our custom resampling really makes a difference."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": false
   },
   "source": [
    "## 1.3 Feature Importance\n",
    "\n",
    "It is worthwile to try and understand what are the most important features according to a random forest classifier. One may even consider to train the final model using only some of the top features."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def feature_importance(n_trees, data):\n",
    "    rf = RandomForestClassifier(n_trees)\n",
    "    rf.fit(data.ix[:,:-1], data.ix[:,-1])\n",
    "    fimp = rf.feature_importances_\n",
    "    results = {}\n",
    "    for idx, name in enumerate(data.ix[:,:-1].columns):\n",
    "        results[name] = fimp[idx]\n",
    "    return results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# We can use an initial random forest to take a look at important features\n",
    "if False: # Set to 'False' to skip this step as it takes some time.\n",
    "    threshold = 0.01\n",
    "    rank = feature_importance(300, data)\n",
    "    print(\"The rank of the features is as follows:\")\n",
    "    flag = 1\n",
    "    count = 0\n",
    "    for a in sorted(rank.keys(), key=rank.get)[::-1]:\n",
    "        aux = rank[a]\n",
    "        count+= 1\n",
    "        if aux < threshold and flag:\n",
    "            flag = 0\n",
    "            print('---'*5)\n",
    "            print('There are {} features above {} '.format(count, threshold))\n",
    "            print('---'*5)\n",
    "        print(a, '-->', rank[a])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Important Remark**: This backs up our data analysis. Notice how all the new features that we created are on the top of the importance rank."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1.4 Ensemble of Random Forests\n",
    "\n",
    "In the final part of the model, we train N_forest-many random forest with n_trees-many trees each. Ultimately, the prediction of the probability of being unsatisfied (TARGET = 0) is made as the geometric mean of the individual predictions.\n",
    "\n",
    "The basic parameters in the model are:\n",
    "- N_forest: number of random forests to train.\n",
    "- n_trees: number of trees for each forest.\n",
    "- a: the fraction of the TARGET == 1 class to be used in the training of each forest.\n",
    "- w: weight of the TARGET == 0 class in the training set."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Cross Validation\n",
    "\n",
    "We perform a grid-search to find the best parameters, we obtain a=0.25, w=1, N_forest=90, n_trees=200."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "cv_data = pd.read_csv(\"cross-validation/cross-val.csv\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>a</th>\n",
       "      <th>w</th>\n",
       "      <th>N_forest</th>\n",
       "      <th>n_trees</th>\n",
       "      <th>auc_roc</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>count</th>\n",
       "      <td>332.000000</td>\n",
       "      <td>332.000000</td>\n",
       "      <td>332.000000</td>\n",
       "      <td>332.000000</td>\n",
       "      <td>332.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>mean</th>\n",
       "      <td>0.717206</td>\n",
       "      <td>1.406627</td>\n",
       "      <td>40.451807</td>\n",
       "      <td>232.078313</td>\n",
       "      <td>0.805536</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>std</th>\n",
       "      <td>0.554718</td>\n",
       "      <td>0.716878</td>\n",
       "      <td>26.442220</td>\n",
       "      <td>136.373686</td>\n",
       "      <td>0.031897</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>min</th>\n",
       "      <td>0.111111</td>\n",
       "      <td>1.000000</td>\n",
       "      <td>10.000000</td>\n",
       "      <td>50.000000</td>\n",
       "      <td>0.720001</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25%</th>\n",
       "      <td>0.250000</td>\n",
       "      <td>1.000000</td>\n",
       "      <td>20.000000</td>\n",
       "      <td>100.000000</td>\n",
       "      <td>0.786379</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>50%</th>\n",
       "      <td>0.500000</td>\n",
       "      <td>1.000000</td>\n",
       "      <td>40.000000</td>\n",
       "      <td>205.000000</td>\n",
       "      <td>0.820972</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>75%</th>\n",
       "      <td>1.000000</td>\n",
       "      <td>2.000000</td>\n",
       "      <td>50.000000</td>\n",
       "      <td>350.000000</td>\n",
       "      <td>0.828311</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>max</th>\n",
       "      <td>5.000000</td>\n",
       "      <td>3.000000</td>\n",
       "      <td>180.000000</td>\n",
       "      <td>500.000000</td>\n",
       "      <td>0.835501</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                a           w    N_forest     n_trees     auc_roc\n",
       "count  332.000000  332.000000  332.000000  332.000000  332.000000\n",
       "mean     0.717206    1.406627   40.451807  232.078313    0.805536\n",
       "std      0.554718    0.716878   26.442220  136.373686    0.031897\n",
       "min      0.111111    1.000000   10.000000   50.000000    0.720001\n",
       "25%      0.250000    1.000000   20.000000  100.000000    0.786379\n",
       "50%      0.500000    1.000000   40.000000  205.000000    0.820972\n",
       "75%      1.000000    2.000000   50.000000  350.000000    0.828311\n",
       "max      5.000000    3.000000  180.000000  500.000000    0.835501"
      ]
     },
     "execution_count": 40,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "cv_data.describe()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>a</th>\n",
       "      <th>w</th>\n",
       "      <th>N_forest</th>\n",
       "      <th>n_trees</th>\n",
       "      <th>auc_roc</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0.25</td>\n",
       "      <td>1</td>\n",
       "      <td>30</td>\n",
       "      <td>200</td>\n",
       "      <td>0.835188</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>0.25</td>\n",
       "      <td>1</td>\n",
       "      <td>60</td>\n",
       "      <td>300</td>\n",
       "      <td>0.835154</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>0.25</td>\n",
       "      <td>1</td>\n",
       "      <td>60</td>\n",
       "      <td>400</td>\n",
       "      <td>0.835114</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>0.25</td>\n",
       "      <td>1</td>\n",
       "      <td>60</td>\n",
       "      <td>500</td>\n",
       "      <td>0.835370</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>0.25</td>\n",
       "      <td>1</td>\n",
       "      <td>90</td>\n",
       "      <td>200</td>\n",
       "      <td>0.835501</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>331</th>\n",
       "      <td>0.25</td>\n",
       "      <td>1</td>\n",
       "      <td>60</td>\n",
       "      <td>150</td>\n",
       "      <td>0.835243</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "        a  w  N_forest  n_trees   auc_roc\n",
       "0    0.25  1        30      200  0.835188\n",
       "5    0.25  1        60      300  0.835154\n",
       "6    0.25  1        60      400  0.835114\n",
       "7    0.25  1        60      500  0.835370\n",
       "8    0.25  1        90      200  0.835501\n",
       "331  0.25  1        60      150  0.835243"
      ]
     },
     "execution_count": 41,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "cv_data[cv_data['auc_roc'] > 0.835]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We compromise with N_forest = 30 and n_trees = 200 for the sake of\n",
    "computer performance."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Model Training"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "a = 0.25 # It can also be >1 to oversample the minority class (1 ie. unhappy)\n",
    "w = 1  # Weight of the majority class in the final sample.\n",
    "N_forest = 90  # Number of forests to train\n",
    "n_trees = 200  # Number of trees in each forest"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "rfs = trainForests(train, a, w, N_forest, n_trees)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# mean_ensemble takes the predicted probabilities of each random forest and computes the geometric mean. See train_forests.py for more details\n",
    "Y_tot = mean_ensemble(rfs, X_test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.8352569907\n"
     ]
    }
   ],
   "source": [
    "temp = roc_auc_score(Y_test,Y_tot['geometric'])\n",
    "rf_classifiers = [rfs, temp]\n",
    "print(temp)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "pickle.dump(rf_classifiers, open(\"models/rf_classifier_param1.dat\", \"wb\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For comparison, we print the score of the prediction on the training set."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Y_prob = mean_ensemble(rfs, X_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.8656689630768486"
      ]
     },
     "execution_count": 48,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "roc_auc_score(Y_train, Y_prob[\"geometric\"])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This indicates that some overfitting might be going on. We are going to try and fix that."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Help with overfitting: 'min_samples_leaf'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "There is another parameter that we can use to prevent overfitting, 'min_samples_leaf'. This is the number of samples that we require each leaf to have in each decision tree of the random forest. The higher the less prone to over fitting. We run a token experiment to see how scores change."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "msl=1 --> 0.7655553447544401\n",
      "msl=2 --> 0.8193590251613391\n",
      "msl=3 --> 0.8196266209379982\n",
      "msl=4 --> 0.820202612100453\n",
      "msl=5 --> 0.8199581530316028\n",
      "msl=6 --> 0.8177373676679066\n",
      "msl=7 --> 0.8187904264591842\n",
      "msl=8 --> 0.8175645195312751\n",
      "msl=9 --> 0.8156901640437729\n",
      "msl=10 --> 0.8177760793302409\n"
     ]
    }
   ],
   "source": [
    "for msl in range(1, 11):\n",
    "    rf = RandomForestClassifier(n_estimators = 300, min_samples_leaf = msl)\n",
    "    rf.fit(train.ix[:,:-1],train['TARGET'])\n",
    "    temp = pd.DataFrame(rf.predict_proba(X_test))[1]\n",
    "    print(\"msl={} --> {}\".format(msl,roc_auc_score(Y_test, temp)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A good choice for the parameter is msl = 5, let's see whether it improves our model."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 52,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "a = 0.25 # It can also be >1 to oversample the minority class (1 ie. unhappy)\n",
    "w = 1  # Weight of the majority class in the final sample.\n",
    "N_forest = 60  # Number of forests to train\n",
    "n_trees = 300  # Number of trees in each forest"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 53,
   "metadata": {
    "collapsed": false,
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "msl=1 --> Test error = 0.8352349261809229, Training error = 0.8694161891445628\n",
      "msl=2 --> Test error = 0.8180729063621319, Training error = 0.8466065649341591\n",
      "msl=3 --> Test error = 0.8121969726125712, Training error = 0.8369061537207224\n",
      "msl=4 --> Test error = 0.8065117756249479, Training error = 0.8297232318854276\n",
      "msl=5 --> Test error = 0.8042192100360492, Training error = 0.8259110773803449\n"
     ]
    }
   ],
   "source": [
    "for msl in range(1,6):\n",
    "    rfs = trainForests(train, a, w, N_forest, n_trees, msl)\n",
    "    Y_tot = mean_ensemble(rfs, X_test)\n",
    "    temp = roc_auc_score(Y_test,Y_tot['geometric'])\n",
    "    rf_classifiers.append([rfs, temp])\n",
    "    Y_prob = mean_ensemble(rfs, X_train)\n",
    "    print(\"msl={} --> Test error = {}, Training error = {}\".format(msl, temp, roc_auc_score(Y_train, Y_prob[\"geometric\"])))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "# classifier = 0 --> 0.790036013132\n",
      "# classifier = 1 --> 0.79141394516\n",
      "# classifier = 2 --> 0.794230923982\n",
      "# classifier = 3 --> 0.794427247412\n",
      "# classifier = 4 --> 0.793491226504\n",
      "# classifier = 5 --> 0.797241178961\n",
      "# classifier = 6 --> 0.790925365608\n",
      "# classifier = 7 --> 0.795520090055\n",
      "# classifier = 8 --> 0.803636672878\n",
      "# classifier = 9 --> 0.80513739875\n",
      "# classifier = 10 --> 0.786201075597\n",
      "# classifier = 11 --> 0.804799094937\n",
      "# classifier = 12 --> 0.793091525769\n",
      "# classifier = 13 --> 0.794867239878\n",
      "# classifier = 14 --> 0.803189231522\n",
      "# classifier = 15 --> 0.8082328081\n",
      "# classifier = 16 --> 0.8053983921\n",
      "# classifier = 17 --> 0.793660293764\n",
      "# classifier = 18 --> 0.800020856896\n",
      "# classifier = 19 --> 0.794469457796\n",
      "# classifier = 20 --> 0.792651702596\n",
      "# classifier = 21 --> 0.794602239927\n",
      "# classifier = 22 --> 0.786603710744\n",
      "# classifier = 23 --> 0.805113302715\n",
      "# classifier = 24 --> 0.799529004346\n",
      "# classifier = 25 --> 0.800931596718\n",
      "# classifier = 26 --> 0.790401742409\n",
      "# classifier = 27 --> 0.806467985173\n",
      "# classifier = 28 --> 0.791616588861\n",
      "# classifier = 29 --> 0.794762503952\n",
      "# classifier = 30 --> 0.793475595118\n",
      "# classifier = 31 --> 0.797255343141\n",
      "# classifier = 32 --> 0.800343698616\n",
      "# classifier = 33 --> 0.801735738392\n",
      "# classifier = 34 --> 0.800212948001\n",
      "# classifier = 35 --> 0.796697748482\n",
      "# classifier = 36 --> 0.79527980688\n",
      "# classifier = 37 --> 0.796106633813\n",
      "# classifier = 38 --> 0.798955778302\n",
      "# classifier = 39 --> 0.798912439299\n",
      "# classifier = 40 --> 0.802046052432\n",
      "# classifier = 41 --> 0.792733978987\n",
      "# classifier = 42 --> 0.797226337609\n",
      "# classifier = 43 --> 0.800307977796\n",
      "# classifier = 44 --> 0.805270180881\n",
      "# classifier = 45 --> 0.795415128405\n",
      "# classifier = 46 --> 0.792893960857\n",
      "# classifier = 47 --> 0.789804533192\n",
      "# classifier = 48 --> 0.794138546444\n",
      "# classifier = 49 --> 0.795392951024\n",
      "# classifier = 50 --> 0.798992402018\n",
      "# classifier = 51 --> 0.796847290618\n",
      "# classifier = 52 --> 0.800955410598\n",
      "# classifier = 53 --> 0.808154481879\n",
      "# classifier = 54 --> 0.802792183043\n",
      "# classifier = 55 --> 0.797059189003\n",
      "# classifier = 56 --> 0.793007725742\n",
      "# classifier = 57 --> 0.789574576889\n",
      "# classifier = 58 --> 0.808024521299\n",
      "# classifier = 59 --> 0.804616766393\n"
     ]
    }
   ],
   "source": [
    "for idx in range(Y_tot.shape[1]-2):\n",
    "    printroc_auc_scorelassifier = {}'.format(idx),'-->',roc_auc_score(Y_test, Y_tot.ix[:,idx]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It doesn't help, at all. We should stick with the default msl=1 value."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "pickle.dump(rf_classifiers, open(\"models/rf_msl_saldo_classifier_param1.dat\", \"wb\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Meta Ensemble"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we are going to use the random forest classifiers trained above to create a meta predictor."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "train = train.reset_index(drop=True)\n",
    "test = test.reset_index(drop=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "train.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X_train, Y_train = train.ix[:,:-1], train['TARGET']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X_train.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "for rf in rfs:\n",
    "    temp = rf.predict(train.ix[:,:-1])\n",
    "    temp = pd.DataFrame(temp)\n",
    "    X_train = pd.concat([X_train, temp], axis=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "train_meta = pd.concat([X_train, train['TARGET']], axis=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "train_meta.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "X_test, Y_test = test.ix[:,:-1], test['TARGET']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "for rf in rfs:\n",
    "    temp = rf.predict(test.ix[:,:-1])\n",
    "    temp = pd.DataFrame(temp)\n",
    "    X_test = pd.concat([X_test, temp], axis=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Y_test = test['TARGET']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X_test.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from sklearn.linear_model import LogisticRegression"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# random regularization\n",
    "cv = {}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "for _ in range(20):\n",
    "    n_trees = np.random.randint(50,300)\n",
    "    lg = RandomForestClassifier(n_estimators = n_trees)\n",
    "    lg.fit(X_train, Y_train)\n",
    "    cv[n_trees] = roc_auc_score(Y_test, lg.predict(X_test))\n",
    "    print(cv[n_trees])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "for key in sorted(cv.keys(), key=cv.get)[::-1]:\n",
    "    print(key,cv[key])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Score analysis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X_test, Y_test = test.ix[:,:-1], test.ix[:,-1]\n",
    "n = 100\n",
    "a = 0.25\n",
    "w = 1\n",
    "N_forest = 5\n",
    "n_trees = 5"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "scores = []\n",
    "for _ in range(n):\n",
    "    rfs = trainForests(train, a, w, N_forest,n_trees)\n",
    "    Y_prob = mean_ensemble(rfs, X_test)\n",
    "    scores.append(roc_auc_score(Y_test,Y_prob))\n",
    "scores = pd.DataFrame(scores)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "scores.describe()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "plt.title(\"Distribution of scores\")\n",
    "plt.hist(scores)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# If desired, transform probabilities into class labels.\n",
    "def threshold(Y_prob, threshold = 0.5):\n",
    "    result = []\n",
    "    for y in Y_prob:\n",
    "        if y <= threshold:\n",
    "            result.append(0)\n",
    "        else:\n",
    "            result.append(1)\n",
    "    return result"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Evaluate class labels\n",
    "Y_pred = threshold(Y_prob, threshold = 0.5)\n",
    "_ = eval_classification(test['TARGET'],Y_pred, print_results = True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Plot feature importance\n",
    "def plot_features(forest):  \n",
    "    importances = forest.feature_importances_\n",
    "    std = np.std([tree.feature_importances_ for tree in forest.estimators_], axis=0)\n",
    "    indices = np.argsort(importances)[::-1]\n",
    "    n=len(indices)\n",
    "    # Plot the feature importances of the forest\n",
    "    plt.figure()\n",
    "    plt.title(\"Feature importances\")\n",
    "    plt.bar(range(n), importances[indices],\n",
    "           color=\"r\", yerr=std[indices], align=\"center\")\n",
    "    plt.xticks(range(n), indices)\n",
    "    plt.xlim([-1, n])\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "train.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Create Submission"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Retrain forest on the whole 'train.csv' data\n",
    "rfs = trainForests(data, a, w, N_forest, n_trees, msl)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "test = pd.read_csv('data/test.csv')\n",
    "test_id = test.ix[:,'ID'].values\n",
    "process(test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "Y_prob = mean_ensemble(rfs,test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "create_submission(test_id, Y_prob['geometric'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Ensemble RF and XGBOOST"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Y_boost = pd.read_csv('../Kaggle_Santander-master/simplexgbtest.csv')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "Y_boost.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Y_rf = pd.read_csv('submissions/rforest_ensemble2.csv')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "Y_rf.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Y_prob = pd.concat([Y_boost,Y_rf.ix[:,'TARGET']], axis=1, ignore_index=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "Y_prob.rename(columns ={0:'ID', 1:'xgb', 2: 'rfe' }, inplace = True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# geometric mean ensemble\n",
    "l = 2 #number of predictors to ensemble\n",
    "temp = Y_prob.ix[:,1:].product(axis=1)\n",
    "temp = temp.apply(lambda x: np.power(x, 1./l))\n",
    "Y_prob['geometric'] = temp"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# arithmetic mean ensemble\n",
    "l = 2 #number of predictors to ensemble\n",
    "temp = Y_prob[['xgb', 'rfe']].mean(axis=1)\n",
    "temp = temp.apply(lambda x: np.power(x, 1./l))\n",
    "Y_prob['arithmetic'] = temp"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# difference column\n",
    "temp = Y_prob['xgb'] - Y_prob['rfe']\n",
    "Y_prob['xgb - rfe'] = temp"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# difference column\n",
    "temp = Y_prob['geometric'] - Y_prob['arithmetic']\n",
    "Y_prob['geo - ari'] = temp"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "Y_prob.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "plt.title('Differences between XGB and RFE')\n",
    "plt.hist(Y_prob['xgb - rfe'])\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "plt.title('Differences between ensembles')\n",
    "plt.hist(Y_prob['geo - ari'])\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "create_submission(test_id, Y_prob['arithmetic'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "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.5.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}