3_Random_Forest_Clustering.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# The Random Forest Algorithm\n",
    "\n",
    "The random forest algorithm is a machine learning method that joins multiple decision tree instances to perform classification or regression. Such combinations are known as ensemble learning models. The term 'random forest' is derived from the inverted forest-like appearance of a multiple decision trees joined at the head.\n",
    "\n",
    "Here we will demonstrate an implementation of both regression (making a prediction based on continuuous data - essentially fancy fitting) and classification (predicting labels of objects based on past experience). Both methods are types of supervised learning, where we train the algorithm on example data.\n",
    "\n",
    "First we need to understand the basics - the decision tree."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Decision Trees\n",
    "\n",
    "Decision Trees (DTs) are a non-parametric supervised learning method used for classification and regression. The goal is to create a model that predicts the value of a target variable by learning simple decision rules inferred from the data features.\n",
    "\n",
    "We can think of a DT as a series of sequential yes/no questions like a flowchart. The final result of the tree is a class label (classification) or value (regression).\n",
    "\n",
    "<img src=\"Images/simple_dtree.png\">\n",
    "\n",
    "A decision tree model is trained and built by determing the questions that lead to the 'best result'. In pracise the best result is determined by the greatest reduction in the *Gini Impurity*. This is simply a measurement of the likelihood of an incorrect classification of a new instance of a random variable - 0 is perfect classification."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1.1 A toy example\n",
    "\n",
    "We start with a very simple example - binary classification in 2D. This is where we can plot two classes of object on a 2D scatter plot, where the x- and y-axis are the features that describe the data. Traditionally in the ML world denote the features as an array *X* and the as a vector *y*.\n",
    "\n",
    "### Create and plot simple data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "\n",
    "# 2D Array of features [x1, x2]\n",
    "X = np.array([[2, 2], \n",
    "              [2, 1],\n",
    "              [2, 3], \n",
    "              [1, 2], \n",
    "              [1, 1],\n",
    "              [3, 3]])\n",
    "\n",
    "# Labels each point as belonging to one of two classes\n",
    "y = np.array([0, 1, 1, 1, 0, 1])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It is easier to understand the data if we visualise it in 2D:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "\n",
    "# Plot formatting\n",
    "plt.rcParams['font.size'] = 18\n",
    "plt.figure(figsize = (8, 8))\n",
    "\n",
    "# Plot each point as the label\n",
    "for x1, x2, label in zip(X[:, 0], X[:, 1], y):\n",
    "    plt.text(x1, x2, str(label), fontsize = 40, color = 'g',\n",
    "             ha='center', va='center')\n",
    "    \n",
    "# Plot formatting\n",
    "plt.grid(None);\n",
    "plt.xlim((0, 3.5));\n",
    "plt.ylim((0, 3.5));\n",
    "plt.xlabel('x1', size = 20); plt.ylabel('x2', size = 20); plt.title('Data', size = 24)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can see from this plot the the two classes are not separable by a single straight line. This means that a simple linear classifier will not work. A decision tree works because it draws many boundaries between the points."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will use a single decision tree to learn how to separate the points. As the tree builds, the algorithm makes branches by maximising the *reduction* in Gini impurity.\n",
    "\n",
    "The default parameters of the scikit-learn decision tree allow it to grow to an arbritrary depth in order to correctly classify all the data. This is essentially menorising the data and leads to overfitting. In real problems we would limit the depth.\n",
    "\n",
    "### Run the decision tree algorithm"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Import the decision tree algorithm from scikit-learn"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Print some stats\n",
    "print(f'Decision tree has {tree.tree_.node_count} nodes with maximum depth {tree.tree_.max_depth}.')\n",
    "print(f'Model Accuracy: {tree.score(X, y)}')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Visualize the decision tree\n",
    "\n",
    "We can plot the decision tree as a flowchart to show each node in the tree. Because the tree is relatively small, we can understand the entire image.\n",
    "\n",
    "First we export the tree as a dot file making sure to label the features and the classes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.tree import export_graphviz\n",
    "\n",
    "# Plot the tree\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A decision tree is an intuitive model: it makes decisions much as we might when faced with a problem by constructing a flowchart of questions. For each of the nodes (except the leaf nodes), the five rows represent:\n",
    "\n",
    " * The node question\n",
    "\n",
    " * gini: the Gini Impurity of the node. The average (weighted by samples) gini impurity decreases with each level of the tree.\n",
    " * samples: number of training observations in the node\n",
    " * value: [number of samples in the first class, number of samples in the second class]\n",
    " * class: the predicted class at this levelp\n",
    " "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercise: Run again with limited depth\n",
    "\n",
    "In practice, we usually want to limit the maximum depth of the decision tree (even in a random forest) so the tree can generalize better to testing data. Although this will lead to reduced accuracy on the training data, it can improve performance on the testing data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Limit maximum depth and train again\n",
    "\n",
    "# Assess the accuracy of the model\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Visualise the tree again\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Our model no longer gets perfect accuracy on the training data. However, it probably would do better on the testing data since we have limited the maximum depth to prevent overfitting.\n",
    "\n",
    "Limiting the depth of a single decision tree is one way we can try to make a less biased model. Another option is to use an entire forest of trees, training each one on a random subsample of the training data. The final model then takes an average of all the individual decision trees to arrive at a classification. This is the idea behind the random forest.\n",
    "<img src=\"Images/Random_forest_diagram_complete.png\">\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1.2 A real example - predicting health\n",
    "\n",
    "The following data set is from the Centers for Disease Control and Prevention (CDC) and includes socioeconomic and lifestyle indicators for hundreds of thousands of individuals. The objective is to predict the overall health of an individual: either 0 for poor health or 1 for good health. We'll limit the data to 100,000 individuals to speed up training.\n",
    "\n",
    "The problem is imbalanced (far more of one label than another) so for assessing performance, we'll use recall, precision, receiver operating characteristic area under the curve (ROC AUC), and also plot the ROC curve. Accuracy is not a useful metric when dealing with an imbalanced problem.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Load the data and clean\n",
    "\n",
    "First step is to load the data into a Pandas *dataframe* and clean any errant values. (A dataframe is a container for tabular data similar to a spreadsheet)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "\n",
    "# Read the data into a dataframe\n",
    "df = pd.read_csv('DATA/2015_health_50k.csv')\n",
    "print(\"Dimensions of data [rows, columns]:\", df.shape)\n",
    "df.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The column '_RFHLTH' gives the ground-truth for health. However, we want to replace values of 2 or above with 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Remove columns with missing values\n",
    "df = df.drop(columns = ['POORHLTH', 'PHYSHLTH', 'GENHLTH', 'PAINACT2', \n",
    "                        'QLMENTL2', 'QLSTRES2', 'QLHLTH2', 'HLTHPLN1', 'MENTHLTH'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Split into train and test set\n",
    "\n",
    "In supervised learning we want to set aside some data as an 'independent' dataset to test our predictions. This is usually some randomly-chosen fraction of the input data. Scikit-learn has a convenience function to make this split easy:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.model_selection import train_test_split\n",
    "\n",
    "# Extract the labels\n",
    "labels = np.array(df.pop('label'))\n",
    "\n",
    "# 30% examples in test data\n",
    "train, test, train_labels, test_labels = train_test_split(df, labels, \n",
    "                                                          stratify = labels,\n",
    "                                                          test_size = 0.3, \n",
    "                                                          random_state = 42)\n",
    "print(\"Train data shape:\", train.shape)\n",
    "print(\"Test data shape:\", test.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Fill in missing values with mean of test \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Grab the list of column names for later use\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Decision Tree on Real Data\n",
    "\n",
    "First, we'll train the decision tree on the data. Let's leave the depth unlimited and see if we get overfitting!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Train tree\n",
    "tree.fit(train, train_labels)\n",
    "print(f'Decision tree has {tree.tree_.node_count} nodes with maximum depth {tree.tree_.max_depth}.')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Assess Decision Tree Performance\n",
    "\n",
    "Given the number of nodes in our decision tree and the maximum depth, we expect it has overfit to the training data. This means it will do much better on the training data than on the testing data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Make probability predictions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.metrics import precision_score, recall_score, roc_auc_score, roc_curve\n",
    "\n",
    "print(f'Train ROC AUC Score: {roc_auc_score(train_labels, train_probs)}')\n",
    "print(f'Test ROC AUC  Score: {roc_auc_score(test_labels, probs)}')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(f'Baseline ROC AUC: {roc_auc_score(test_labels, [1 for _ in range(len(test_labels))])}')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Evaluate the Decision Tree\n",
    "\n",
    "We'll write a short function that calculates a number of metrics for the baseline (guessing the most common label in the training data), the testing predictions, and the training predictions. The function also plots the ROC curve where a better model is to the left and towards the top."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def evaluate_model(predictions, probs, train_predictions, train_probs):\n",
    "    \"\"\"Compare machine learning model to baseline performance.\n",
    "    Computes statistics and shows ROC curve.\"\"\"\n",
    "    \n",
    "    baseline = {}\n",
    "    \n",
    "    baseline['recall'] = recall_score(test_labels, [1 for _ in range(len(test_labels))])\n",
    "    baseline['precision'] = precision_score(test_labels, [1 for _ in range(len(test_labels))])\n",
    "    baseline['roc'] = 0.5\n",
    "    \n",
    "    results = {}\n",
    "    \n",
    "    results['recall'] = recall_score(test_labels, predictions)\n",
    "    results['precision'] = precision_score(test_labels, predictions)\n",
    "    results['roc'] = roc_auc_score(test_labels, probs)\n",
    "    \n",
    "    train_results = {}\n",
    "    train_results['recall'] = recall_score(train_labels, train_predictions)\n",
    "    train_results['precision'] = precision_score(train_labels, train_predictions)\n",
    "    train_results['roc'] = roc_auc_score(train_labels, train_probs)\n",
    "    \n",
    "    for metric in ['recall', 'precision', 'roc']:\n",
    "        print(f'{metric.capitalize()} Baseline: {round(baseline[metric], 2)} Test: {round(results[metric], 2)} Train: {round(train_results[metric], 2)}')\n",
    "    \n",
    "    # Calculate false positive rates and true positive rates\n",
    "    base_fpr, base_tpr, _ = roc_curve(test_labels, [1 for _ in range(len(test_labels))])\n",
    "    model_fpr, model_tpr, _ = roc_curve(test_labels, probs)\n",
    "\n",
    "    plt.figure(figsize = (8, 6))\n",
    "    plt.rcParams['font.size'] = 16\n",
    "    \n",
    "    # Plot both curves\n",
    "    plt.plot(base_fpr, base_tpr, 'b', label = 'baseline')\n",
    "    plt.plot(model_fpr, model_tpr, 'r', label = 'model')\n",
    "    plt.legend();\n",
    "    plt.xlabel('False Positive Rate'); plt.ylabel('True Positive Rate'); plt.title('ROC Curves');"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from collections import Counter\n",
    "print(Counter(probs))\n",
    "print(Counter(predictions))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "evaluate_model(predictions, probs, train_predictions, train_probs)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "There we can see the problem with a single decision tree where the maximum depth is not limited: severe overfitting to the training data.\n",
    "\n",
    "Another method to inspect the performance of a classification model is by making a confusion matrix.\n",
    "\n",
    "### Confusion Matrix"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.metrics import confusion_matrix\n",
    "import itertools\n",
    "\n",
    "def plot_confusion_matrix(cm, classes,\n",
    "                          normalize=False,\n",
    "                          title='Confusion matrix',\n",
    "                          cmap=plt.cm.Oranges):\n",
    "    \"\"\"\n",
    "    This function prints and plots the confusion matrix.\n",
    "    Normalization can be applied by setting `normalize=True`.\n",
    "    Source: http://scikit-learn.org/stable/auto_examples/model_selection/plot_confusion_matrix.html\n",
    "    \"\"\"\n",
    "    if normalize:\n",
    "        cm = cm.astype('float') / cm.sum(axis=1)[:, np.newaxis]\n",
    "        print(\"Normalized confusion matrix\")\n",
    "    else:\n",
    "        print('Confusion matrix, without normalization')\n",
    "\n",
    "    print(cm)\n",
    "\n",
    "    plt.figure(figsize = (10, 10))\n",
    "    plt.imshow(cm, interpolation='nearest', cmap=cmap)\n",
    "    plt.title(title, size = 24)\n",
    "    plt.colorbar(aspect=4)\n",
    "    tick_marks = np.arange(len(classes))\n",
    "    plt.xticks(tick_marks, classes, rotation=45, size = 14)\n",
    "    plt.yticks(tick_marks, classes, size = 14)\n",
    "\n",
    "    fmt = '.2f' if normalize else 'd'\n",
    "    thresh = cm.max() / 2.\n",
    "    \n",
    "    # Labeling the plot\n",
    "    for i, j in itertools.product(range(cm.shape[0]), range(cm.shape[1])):\n",
    "        plt.text(j, i, format(cm[i, j], fmt), fontsize = 20,\n",
    "                 horizontalalignment=\"center\",\n",
    "                 color=\"white\" if cm[i, j] > thresh else \"black\")\n",
    "        \n",
    "#    plt.grid(None)\n",
    "#    plt.tight_layout()\n",
    "    plt.ylabel('True label', size = 18)\n",
    "    plt.xlabel('Predicted label', size = 18)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cm = confusion_matrix(test_labels, predictions)\n",
    "plot_confusion_matrix(cm, classes = ['Poor Health', 'Good Health'],\n",
    "                      title = 'Health Confusion Matrix')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This shows the classifications predicted by the model on the test data along with the real labels. We can see that our model has many false negatives (predicted good health but actually poor health) and false positives (predicted poor health but actually good health)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Feature Importances\n",
    "Finally, we can take a look at the features considered most important by the Decision Tree. The values are computed by summing the reduction in Gini Impurity over all of the nodes of the tree in which the feature is used."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fi = pd.DataFrame({'feature': features,\n",
    "                   'importance': tree.feature_importances_}).\\\n",
    "                    sort_values('importance', ascending = False)\n",
    "fi.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualize Full Tree\n",
    "\n",
    "As before, we can look at the decision tree on the data. This time, we have to limit the maximum depth otherwise the tree will be too large and cannot be converted and displayed as an image."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Save tree as dot file\n",
    "export_graphviz(tree, 'TMP/tree_real_data.dot', rounded = True, \n",
    "                feature_names = features, max_depth = 6,\n",
    "                class_names = ['poor health', 'good health'], filled = True)\n",
    "\n",
    "# Convert to png\n",
    "!dot -Tpng TMP/tree_real_data.dot -o TMP/tree_real_data.png -Gdpi=200\n",
    "\n",
    "# Visualize\n",
    "Image(filename='TMP/tree_real_data.png')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can see that our model is extremely deep and has many nodes. To reduce the variance of our model, we could limit the maximum depth or the number of leaf nodes. Another method to reduce the variance is to use more trees, each one trained on a random sampling of the observations. This is where the random forest comes into play."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. Random Forest\n",
    "\n",
    "This takes the idea of a single decision tree, and creates an ensemble model out of hundreds or thousands of trees to reduce the variance. Each tree is trained on a random set of the observations, and for each split of a node, only a subset of the features are used for making a split. When making predictions, the random forest averages the predictions for each of the individual decision trees for each data point in order to arrive at a final classification.\n",
    "\n",
    "### 2.1 Random Forest with fixed parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.ensemble import RandomForestClassifier\n",
    "\n",
    "# Create the model with 100 trees\n",
    "model = RandomForestClassifier(n_estimators=100, \n",
    "                               random_state=42, \n",
    "                               max_features = 'sqrt',\n",
    "                               n_jobs=-1, verbose = 1)\n",
    "\n",
    "# Fit on training data\n",
    "model.fit(train, train_labels)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can see how many nodes there are for each tree on average and the maximum depth of each tree. There were 100 trees in the forest."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n_nodes = []\n",
    "max_depths = []\n",
    "\n",
    "for ind_tree in model.estimators_:\n",
    "    n_nodes.append(ind_tree.tree_.node_count)\n",
    "    max_depths.append(ind_tree.tree_.max_depth)\n",
    "    \n",
    "print(f'Average number of nodes {int(np.mean(n_nodes))}')\n",
    "print(f'Average maximum depth {int(np.mean(max_depths))}')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Results\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "train_rf_predictions = model.predict(train)\n",
    "train_rf_probs = model.predict_proba(train)[:, 1]\n",
    "\n",
    "rf_predictions = model.predict(test)\n",
    "rf_probs = model.predict_proba(test)[:, 1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "evaluate_model(rf_predictions, rf_probs, train_rf_predictions, train_rf_probs)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cm = confusion_matrix(test_labels, rf_predictions)\n",
    "plot_confusion_matrix(cm, classes = ['Poor Health', 'Good Health'],\n",
    "                      title = 'Health Confusion Matrix')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compared to the single decision tree, the model has fewer false postives although more false negatives. Overall, the random forest does significantly better than a single decision tree. This is what we expected!"
   ]
  }
 ],
 "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.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}