pca-svm_classifier_skopt.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# PCA - SVM Classifier\n",
    "This notebook shows how to train a Support Vector Machine classifier on the grayscale Rock-Paper-Scissors images using Principal Component Analysis for dimensionality reduction and Bayesian optimization of the model hyperparameters.\n",
    "\n",
    "The model is trained using 75% of the dataset with a five fold cross-validation and achieves a f1 (micro) score of 97.6% on a test set of 547 images (the remainig 25% of the dataset)\n",
    "\n",
    "**Author**: Julien de la Bruère-Terreault <drgfreeman@tuta.io>  \n",
    "**Source**: https://github.com/DrGFreeman/rps-cv-data-science  \n",
    "**License**: MIT"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Import packages\n",
    "In addition to the common Numpy, and Scikit-Learn packages, this model also requires the [Scikit-Optimize](https://github.com/scikit-optimize/scikit-optimize) package for optimization and the [rpscv](https://github.com/DrGFreeman/rps-cv) package for preprocessing of the images. The instructions to install these packages are given in the README of the [GitHub repository](https://github.com/DrGFreeman/rps-cv-data-science)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "\n",
    "from sklearn.model_selection import train_test_split\n",
    "from sklearn.pipeline import Pipeline\n",
    "from sklearn.decomposition import PCA\n",
    "from sklearn.svm import SVC\n",
    "from sklearn.model_selection import StratifiedKFold\n",
    "from sklearn.metrics import classification_report\n",
    "from sklearn.metrics import confusion_matrix\n",
    "\n",
    "from skopt import BayesSearchCV\n",
    "from skopt.space import Real, Integer\n",
    "\n",
    "from rpscv import imgproc"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Generate feature and labels from images\n",
    "The `generateGrayFeatures` function of the `rpscv.imgproc` module reads the images, removes the green background, converts them to grayscale and rescales their values in the range of 0 (black) to 1 (white). The function outputs an array where each row represents the 60,000 pixel values of an image (200x300 pixels flattened) and a vector of label values in the form of integers (0 = *rock*, 1 = *paper*, 2 = *scissors*)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Completed processing 2188 images\n"
     ]
    }
   ],
   "source": [
    "X, y = imgproc.generateGrayFeatures()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Train-test split\n",
    "We split the dataset into training and testing sets using 75% of the images for training and 25% for testing."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=.25,\n",
    "                                                    stratify=y, random_state=42)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Model definition\n",
    "The model consists of a Principal Component Analysis (PCA) transformer followed by a Support Vector Machine classifier. The PCA transformer will reduce the number of features from 60,000 (the number of pixels in the original images) to a reduced number of values. For more information on the PCA transformation, see the [PCA visualization](pca_visualization.ipynb) notebook in this repository.\n",
    "\n",
    "The SVC classifier will take the transformed features from the PCA, along with the image labels, to fit a classification model and perform the predictions.\n",
    "\n",
    "The PCA transformer and SVC classifier are assembled into a Scikit-Learn `Pipeline` object. This will allow to perform cross-validation on the training data during the model hyperparameter tuning without leaking data from the hold-out fold into the model.\n",
    "\n",
    "The `rbf` kernel is selected for the support vector classifier as it allows more complex, non-linear decision boundaries."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "pipeline = Pipeline([('pca', PCA()),\n",
    "                     ('clf', SVC(kernel='rbf'))])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, we define the range of the model hyperparameters into which to perform the search. Since we will be using the `BayesSearchCV` class of the Scikit-Optimize package, we define the hyperparameter ranges using the `Integer` and `Real` classes of the same package.\n",
    "\n",
    "We define the range of principal components (`n_components`) for the `PCA` transformer as integers between 20 and 100.\n",
    "\n",
    "We define the `gamma` and `C` parameters of the `SVC` as reals with the ranges .0001-.01 and 1-3000 respectively."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "opt_params = {'pca__n_components': Integer(20, 100),\n",
    "              'clf__gamma': Real(.0001, .01, prior='log-uniform'),\n",
    "              'clf__C': Real(1, 3000, prior='log-uniform')}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We define the model as the `BayesSearchCV` cross-validation object, passing the `pipeline` object defined previously as estimator. We pass the hyperparameter dictionary as search space and specify 100 iterations for the hyperparameter search.\n",
    "\n",
    "We use [5 fold cross-validation](https://en.wikipedia.org/wiki/Cross-validation_(statistics)#k-fold_cross-validation) which means that at each search iteration, the training set will be split in five and the model fitted to 4 of the folds (80% of the training data) and evaluated on the remaining 20% of the training data. This will be repeated 5 times, changing the validation set each time, and the average score at each iteration will be kept. The use of `StratifiedKFold` ensures that each class will be equally represented in each of the cross-validation folds.\n",
    "\n",
    "The scoring function used is the [F1 score](https://en.wikipedia.org/wiki/F1_score) calculated globally (micro option). This metric balances the precision and recall taking into account class imbalance if present."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "model = BayesSearchCV(pipeline,\n",
    "                      search_spaces=opt_params,\n",
    "                      n_iter=100,\n",
    "                      cv=StratifiedKFold(n_splits=5),\n",
    "                      scoring='f1_micro',\n",
    "                      n_jobs=-1,\n",
    "                      return_train_score=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Fit model and evaluate on training set\n",
    "With the model defined, we call the `fit` method on the model which performs the hyperparameter optimization."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CPU times: user 6min 22s, sys: 9min 25s, total: 15min 48s\n",
      "Wall time: 35min 25s\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "BayesSearchCV(cv=StratifiedKFold(n_splits=5, random_state=None, shuffle=False),\n",
       "       error_score='raise',\n",
       "       estimator=Pipeline(memory=None,\n",
       "     steps=[('pca', PCA(copy=True, iterated_power='auto', n_components=None, random_state=None,\n",
       "  svd_solver='auto', tol=0.0, whiten=False)), ('clf', SVC(C=1.0, cache_size=200, class_weight=None, coef0=0.0,\n",
       "  decision_function_shape='ovr', degree=3, gamma='auto_deprecated',\n",
       "  kernel='rbf', max_iter=-1, probability=False, random_state=None,\n",
       "  shrinking=True, tol=0.001, verbose=False))]),\n",
       "       fit_params=None, iid=True, n_iter=100, n_jobs=-1, n_points=1,\n",
       "       optimizer_kwargs=None, pre_dispatch='2*n_jobs', random_state=None,\n",
       "       refit=True, return_train_score=True, scoring='f1_micro',\n",
       "       search_spaces={'pca__n_components': Integer(low=20, high=100), 'clf__gamma': Real(low=0.0001, high=0.01, prior='log-uniform', transform='identity'), 'clf__C': Real(low=1, high=3000, prior='log-uniform', transform='identity')},\n",
       "       verbose=0)"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "%%time\n",
    "model.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We print the best set of hyperparamters found during the optimization as well as the model F1 score obtained on the training set with the k-fold cross-validation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'clf__C': 289.6483653005734,\n",
       " 'clf__gamma': 0.0006175507672838874,\n",
       " 'pca__n_components': 67}"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model.best_params_"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.9847653869591713"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model.best_score_"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Model evaluation on test set\n",
    "By default, the `BayesSearchCV` object, as well as the Scikit-Learn `GridSearchCV` of `RandomSearchCV` objects, retrain the model, with the best hyperparameters, on the full training set at the end of the hyperparameter optimization. We therefore don't need to perform this test manually prior evaluating the model on the test set.\n",
    "\n",
    "We use the model `score` method to evaluate the model F1 score on the test set which it has never seen."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.9762340036563071"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model.score(X_test, y_test)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The model achieves a score of 0.976 which is excellent for a \"traditional\" (i.e. non deep-learning) model on such a small dataset with no data augmentation used. There is also only a small loss of performance compared to the model score on the training set which indicates that the model is not overfitting the data.\n",
    "\n",
    "To further evaluate the model performance, we use the model's `predict` method to predict the labels of the test set."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "y_pred = model.predict(X_test)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For convenience, we import the `gestureTxt` dictionary from the `rpscv.utils` module. This dictionary provides the text labels matching the integer values returned by the `generateGrayFeatures` function."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "from rpscv.utils import gestureTxt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We print the classification report using the predicted labels `y_pred` and the true labels `y_test`. This report provides the precision, recall and f1-score values for each class as well as the average values for each score using different weightings."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "              precision    recall  f1-score   support\n",
      "\n",
      "        rock       0.96      0.99      0.98       182\n",
      "       paper       0.97      0.97      0.97       178\n",
      "    scissors       1.00      0.97      0.99       187\n",
      "\n",
      "   micro avg       0.98      0.98      0.98       547\n",
      "   macro avg       0.98      0.98      0.98       547\n",
      "weighted avg       0.98      0.98      0.98       547\n",
      "\n"
     ]
    }
   ],
   "source": [
    "print(classification_report(y_test, y_pred, target_names=gestureTxt.values()))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, we print the confusion matrix of predictions on the test set which helps undertand where the models predictions are wrong.\n",
    "\n",
    "For instance, looking at the columns, we can see that the models predicted *rock* for six of the *paper* images and one *scissors* image, which explains the lower precision score of 96% of the *rock* class reported in the classification report. In opposition, all predicted *scissors* are indeed *scissors*, explaining the 100% precision score on the *scissors* class.\n",
    "\n",
    "Looking at the rows, we can see that out of the 182 *rock* images, only two we incorrectly predicted as *paper*. This shows in the 99% recall score for the *rock* class in the classification report compared to the 97% recall scores of the *paper* and *scissors* classes which have mode incorrectly predicted images."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "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 tr th {\n",
       "        text-align: left;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th colspan=\"3\" halign=\"left\">predicted label</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th>rock</th>\n",
       "      <th>paper</th>\n",
       "      <th>scissors</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th rowspan=\"3\" valign=\"top\">true label</th>\n",
       "      <th>rock</th>\n",
       "      <td>180</td>\n",
       "      <td>2</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>paper</th>\n",
       "      <td>6</td>\n",
       "      <td>172</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>scissors</th>\n",
       "      <td>1</td>\n",
       "      <td>4</td>\n",
       "      <td>182</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                    predicted label               \n",
       "                               rock paper scissors\n",
       "true label rock                 180     2        0\n",
       "           paper                  6   172        0\n",
       "           scissors               1     4      182"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "conf_matrix = pd.DataFrame(confusion_matrix(y_test, y_pred))\n",
    "conf_matrix.index = pd.MultiIndex.from_tuples([('true label', label) for label in gestureTxt.values()])\n",
    "conf_matrix.columns = pd.MultiIndex.from_tuples([('predicted label', label) for label in gestureTxt.values()])\n",
    "conf_matrix"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python [conda env:rps-cv]",
   "language": "python",
   "name": "conda-env-rps-cv-py"
  },
  "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.7.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}