class06_Scikit-Learn.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "# Data Mining:<br>Statistical Modeling and Learning from Data\n",
    "\n",
    "## Dr. Ciro Cattuto<br>Dr. Laetitia Gauvin<br>Dr. André Panisson\n",
    "\n",
    "### Exercises - Scikit-Learn\n",
    "\n",
    "---"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%pylab inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Scikit-Learn: Introduction\n",
    "\n",
    "Scikit-learn is a collection of tools for machine learning written in Python:\n",
    "[http://scikit-learn.org](http://scikit-learn.org)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Representation of Data in Scikit-learn\n",
    "\n",
    "Most machine learning algorithms implemented in scikit-learn expect data to be stored in a\n",
    "**two-dimensional array or matrix**.  The arrays can be\n",
    "either ``numpy`` arrays, or in some cases ``scipy.sparse`` matrices.\n",
    "The size of the array is expected to be `[n_samples, n_features]`\n",
    "\n",
    "- **n_samples:**   The number of samples: each sample is an item to process (e.g. classify).\n",
    "  A sample can be a document, a picture, a sound, a video, an astronomical object,\n",
    "  a row in database or CSV file,\n",
    "  or whatever you can describe with a fixed set of quantitative traits.\n",
    "- **n_features:**  The number of features or distinct traits that can be used to describe each\n",
    "  item in a quantitative manner.  Features are generally real-valued, but may be boolean or\n",
    "  discrete-valued in some cases.\n",
    "\n",
    "The number of features must be fixed in advance. However it can be very high dimensional\n",
    "(e.g. millions of features) with most of them being zeros for a given sample. This is a case\n",
    "where `scipy.sparse` matrices can be useful, in that they are\n",
    "much more memory-efficient than numpy arrays.\n",
    "\n",
    "\n",
    "A classification algorithm, for example, expects the data to be represented as a **feature matrix** and a **label vector**:\n",
    "\n",
    "$$\n",
    "{\\rm feature~matrix:~~~} {\\bf X}~=~\\left[\n",
    "\\begin{matrix}\n",
    "x_{11} & x_{12} & \\cdots & x_{1D}\\\\\n",
    "x_{21} & x_{22} & \\cdots & x_{2D}\\\\\n",
    "x_{31} & x_{32} & \\cdots & x_{3D}\\\\\n",
    "\\vdots & \\vdots & \\ddots & \\vdots\\\\\n",
    "\\vdots & \\vdots & \\ddots & \\vdots\\\\\n",
    "x_{N1} & x_{N2} & \\cdots & x_{ND}\\\\\n",
    "\\end{matrix}\n",
    "\\right]\n",
    "$$\n",
    "\n",
    "$$\n",
    "{\\rm label~vector:~~~} {\\bf y}~=~ [y_1, y_2, y_3, \\cdots y_N]\n",
    "$$\n",
    "\n",
    "Here there are $N$ samples and $D$ features."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## A Simple Example: the Iris Dataset\n",
    "\n",
    "As an example of a simple dataset, we're going to take a look at the\n",
    "iris data stored by scikit-learn.\n",
    "The data consists of measurements of three different species of irises.\n",
    "There are three species of iris in the dataset, which we can picture here:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from IPython.core.display import Image, display\n",
    "display(Image('http://mirlab.org/jang/books/dcpr/image/Iris-setosa-10_1.jpg', width=200, height=200))\n",
    "print (\"Iris Setosa\\n\")\n",
    "\n",
    "display(Image('http://mirlab.org/jang/books/dcpr/image/Iris-versicolor-21_1.jpg', width=200, height=200))\n",
    "print (\"Iris Versicolor\\n\")\n",
    "\n",
    "display(Image('http://mirlab.org/jang/books/dcpr/image/Iris-virginica-3_1.jpg', width=200, height=200))\n",
    "print (\"Iris Virginica\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Loading the Iris Data with Scikit-Learn\n",
    "\n",
    "Scikit-learn has a very straightforward set of data on these iris species.  The data consist of\n",
    "the following:\n",
    "\n",
    "- Features in the Iris dataset:\n",
    "\n",
    "  1. sepal length in cm\n",
    "  2. sepal width in cm\n",
    "  3. petal length in cm\n",
    "  4. petal width in cm\n",
    "\n",
    "- Target classes to predict:\n",
    "\n",
    "  1. Iris Setosa\n",
    "  2. Iris Versicolour\n",
    "  3. Iris Virginica\n",
    "  \n",
    "``scikit-learn`` embeds a copy of the iris CSV file along with a helper function to load it into numpy arrays:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.datasets import load_iris\n",
    "iris = load_iris()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "iris.keys()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n_samples, n_features = iris.data.shape\n",
    "print (n_samples, n_features)\n",
    "print (iris.data[0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print (iris.data.shape)\n",
    "print (iris.target.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print (iris.target)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print (iris.target_names)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This data is four dimensional, but we can visualize two of the dimensions\n",
    "at a time using a simple scatter-plot:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "x_index = 2\n",
    "y_index = 3\n",
    "\n",
    "# this formatter will label the colorbar with the correct target names\n",
    "formatter = plt.FuncFormatter(lambda i, *args: iris.target_names[int(i)])\n",
    "\n",
    "plt.scatter(iris.data[:, x_index], iris.data[:, y_index],\n",
    "            c=iris.target)\n",
    "plt.colorbar(ticks=[0, 1, 2], format=formatter)\n",
    "plt.xlabel(iris.feature_names[x_index])\n",
    "plt.ylabel(iris.feature_names[y_index]);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Other Available Data\n",
    "They come in three flavors:\n",
    "\n",
    "- **Packaged Data:** these small datasets are packaged with the scikit-learn installation,\n",
    "  and can be downloaded using the tools in ``sklearn.datasets.load_*``\n",
    "- **Downloadable Data:** these larger datasets are available for download, and scikit-learn\n",
    "  includes tools which streamline this process.  These tools can be found in\n",
    "  ``sklearn.datasets.fetch_*``\n",
    "- **Generated Data:** there are several datasets which are generated from models based on a\n",
    "  random seed.  These are available in the ``sklearn.datasets.make_*``"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn import datasets"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "datasets.make_classification()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### The Scikit-learn Estimator Object\n",
    "\n",
    "Every algorithm is exposed in scikit-learn via an ''Estimator'' object. For instance a linear regression is implemented as so:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.linear_model import LinearRegression"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Estimator parameters**: All the parameters of an estimator can be set when it is instantiated, and have suitable default values:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model = LinearRegression(normalize=True)\n",
    "print (model.normalize)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print (model)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Estimated Model parameters**: When data is *fit* with an estimator, parameters are estimated from the data at hand. All the estimated parameters are attributes of the estimator object ending by an underscore:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "x = np.array([0, 1, 2])\n",
    "y = np.array([0, 1, 2])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.plot(x, y, 'o')\n",
    "plt.xlim(-0.5, 2.5)\n",
    "plt.ylim(-0.5, 2.5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# The input data for sklearn is 2D: (samples == 3 x features == 1)\n",
    "X = x[:, np.newaxis]\n",
    "print (X)\n",
    "print (y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model.fit(X, y) \n",
    "print (model.coef_)\n",
    "print (model.intercept_)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The model found a line with a slope 1 and intercept 0, as we'd expect."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Supervised vs. Unsupervised Learning\n",
    "\n",
    "The algorithms of machine learning are generally split into two basic categories: **Supervised** and **Unsupervised** learning.\n",
    "\n",
    "### Supervised Learning\n",
    "**Supervised** learning concerns **labeled** data, and the construction of models that can be used to predict labels for new, unlabeled data.\n",
    "\n",
    "*Example:* given a set of *labeled* hand-written digits, create an algorithm which will predict the label of a new instance for which a label is not known.\n",
    "\n",
    "### Unsupervised Learning\n",
    "**Unsupervised** learning concerns **unlabeled** data, and finding structure in the data such as clusters, important dimensions, etc.\n",
    "\n",
    "*Example:* given a set of *unlabeled* digits, determine which digits are related."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Supervised Learning: Classification and Regression\n",
    "\n",
    "In **Supervised Learning**, we have a dataset consisting of both features and labels.\n",
    "The task is to construct an estimator which is able to predict the label of an object\n",
    "given the set of features. A relatively simple example is predicting the species of \n",
    "iris given a set of measurements of its flower. This is a relatively simple task. \n",
    "Some more complicated examples are:\n",
    "\n",
    "- given a multicolor image of an object through a telescope, determine\n",
    "  whether that object is a star, a quasar, or a galaxy.\n",
    "- given a photograph of a person, identify the person in the photo.\n",
    "- given a list of movies a person has watched and their personal rating\n",
    "  of the movie, recommend a list of movies they would like\n",
    "  (So-called *recommender systems*: a famous example is the [Netflix Prize](http://en.wikipedia.org/wiki/Netflix_prize)).\n",
    "\n",
    "What these tasks have in common is that there is one or more unknown\n",
    "quantities associated with the object which needs to be determined from other\n",
    "observed quantities.\n",
    "\n",
    "Supervised learning is further broken down into two categories, **classification** and **regression**.\n",
    "In classification, the label is discrete, while in regression, the label is continuous."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Classification Example\n",
    "K nearest neighbors (kNN) is one of the simplest learning strategies: given a new, unknown observation, look up in your reference database which ones have the closest features and assign the predominant class.\n",
    "\n",
    "Let's try it out on our iris classification problem:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn import neighbors, datasets\n",
    "\n",
    "iris = datasets.load_iris()\n",
    "X, y = iris.data, iris.target\n",
    "\n",
    "# create the model\n",
    "knn = neighbors.KNeighborsClassifier(n_neighbors=1)\n",
    "\n",
    "# fit the model\n",
    "knn.fit(X, y)\n",
    "\n",
    "# What kind of iris has 3cm x 5cm sepal and 4cm x 2cm petal?\n",
    "# call the \"predict\" method:\n",
    "result = knn.predict([[3, 5, 4, 2],])\n",
    "\n",
    "print (iris.target_names[result])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.svm import SVC"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model = SVC()\n",
    "model.fit(X, y)\n",
    "result = model.predict([[3, 5, 4, 2],])\n",
    "print (iris.target_names[result])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Regression Example\n",
    "\n",
    "Simplest possible regression is fitting a line to data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print (X.squeeze())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create some simple data\n",
    "import numpy as np\n",
    "np.random.seed(0)\n",
    "X = np.random.random(size=(20, 1))\n",
    "y = 3 * X.squeeze() + 2 + np.random.normal(size=20)\n",
    "\n",
    "# Fit a linear regression to it\n",
    "from sklearn.linear_model import LinearRegression\n",
    "model = LinearRegression(fit_intercept=True)\n",
    "model.fit(X, y)\n",
    "print (\"Model coefficient: %.5f, and intercept: %.5f\"\n",
    "       % (model.coef_, model.intercept_))\n",
    "\n",
    "# Plot the data and the model prediction\n",
    "X_test = np.linspace(0, 1, 100)[:, np.newaxis]\n",
    "y_test = model.predict(X_test)\n",
    "\n",
    "plt.plot(X.squeeze(), y, 'o')\n",
    "plt.plot(X_test.squeeze(), y_test);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Unsupervised Learning: Dimensionality Reduction and Clustering\n",
    "\n",
    "**Unsupervised Learning** addresses a different sort of problem. Here the data has no labels,\n",
    "and we are interested in finding similarities between the objects in question. In a sense,\n",
    "you can think of unsupervised learning as a means of discovering labels from the data itself.\n",
    "Unsupervised learning comprises tasks such as *dimensionality reduction*, *clustering*, and\n",
    "*density estimation*. For example, in the iris data discussed above, we can used unsupervised\n",
    "methods to determine combinations of the measurements which best display the structure of the\n",
    "data. As we'll see below, such a projection of the data can be used to visualize the\n",
    "four-dimensional dataset in two dimensions. Some more involved unsupervised learning problems are:\n",
    "\n",
    "- given detailed observations of distant galaxies, determine which features or combinations of\n",
    "  features best summarize the information.\n",
    "- given a mixture of two sound sources (for example, a person talking over some music),\n",
    "  separate the two (this is called the [blind source separation](http://en.wikipedia.org/wiki/Blind_signal_separation) problem).\n",
    "- given a video, isolate a moving object and categorize in relation to other moving objects which have been seen.\n",
    "\n",
    "Sometimes the two may even be combined: e.g. Unsupervised learning can be used to find useful\n",
    "features in heterogeneous data, and then these features can be used within a supervised\n",
    "framework."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Dimensionality Reduction: PCA\n",
    "\n",
    "Principle Component Analysis (PCA) is a dimension reduction technique that can find the combinations of variables that explain the most variance.\n",
    "\n",
    "Consider the iris dataset. It cannot be visualized in a single 2D plot, as it has 4 features. We are going to extract 2 combinations of sepal and petal dimensions to visualize it:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X, y = iris.data, iris.target\n",
    "from sklearn.decomposition import PCA\n",
    "pca = PCA(n_components=2)\n",
    "pca.fit(X)\n",
    "X_reduced = pca.transform(X)\n",
    "print (\"Reduced dataset shape:\", X_reduced.shape)\n",
    "\n",
    "import pylab as pl\n",
    "pl.scatter(X_reduced[:, 0], X_reduced[:, 1], c=y)\n",
    "\n",
    "print (\"Meaning of the 2 components:\")\n",
    "for component in pca.components_:\n",
    "    print (\" + \".join(\"%.3f x %s\" % (value, name)\n",
    "                     for value, name in zip(component,\n",
    "                                            iris.feature_names)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Clustering: K-means\n",
    "\n",
    "Clustering groups together observations that are homogeneous with respect to a given criterion, finding ''clusters'' in the data.\n",
    "\n",
    "Note that these clusters will uncover relevent hidden structure of the data only if the criterion used highlights it."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "from sklearn.cluster import KMeans\n",
    "k_means = KMeans(n_clusters=3, random_state=0) # Fixing the RNG in kmeans\n",
    "k_means.fit(X)\n",
    "y_pred = k_means.predict(X)\n",
    "\n",
    "pl.scatter(X_reduced[:, 0], X_reduced[:, 1], c=y_pred);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Scikit-learn's estimator interface\n",
    "\n",
    "Scikit-learn strives to have a uniform interface across all methods. Given a scikit-learn *estimator*\n",
    "object named `model`, the following methods are available:\n",
    "\n",
    "- Available in **all Estimators**\n",
    "  + `model.fit()` : fit training data. For supervised learning applications,\n",
    "    this accepts two arguments: the data `X` and the labels `y` (e.g. `model.fit(X, y)`).\n",
    "    For unsupervised learning applications, this accepts only a single argument,\n",
    "    the data `X` (e.g. `model.fit(X)`).\n",
    "- Available in **supervised estimators**\n",
    "  + `model.predict()` : given a trained model, predict the label of a new set of data.\n",
    "    This method accepts one argument, the new data `X_new` (e.g. `model.predict(X_new)`),\n",
    "    and returns the learned label for each object in the array.\n",
    "  + `model.predict_proba()` : For classification problems, some estimators also provide\n",
    "    this method, which returns the probability that a new observation has each categorical label.\n",
    "    In this case, the label with the highest probability is returned by `model.predict()`.\n",
    "  + `model.score()` : for classification or regression problems, most (all?) estimators implement\n",
    "    a score method.  Scores are between 0 and 1, with a larger score indicating a better fit.\n",
    "- Available in **unsupervised estimators**\n",
    "  + `model.transform()` : given an unsupervised model, transform new data into the new basis.\n",
    "    This also accepts one argument `X_new`, and returns the new representation of the data based\n",
    "    on the unsupervised model.\n",
    "  + `model.fit_transform()` : some estimators implement this method,\n",
    "    which more efficiently performs a fit and a transform on the same input data."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Machine Learning with Scikit-Learn: Validation and Model Selection\n",
    "\n",
    "This section focuses on **validation** and **model selection**."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Supervised Learning Example: Classifying Digits\n",
    "\n",
    "Features can be any **uniformly measured** numerical observation of the data. For example, in the digits data, the features are the brightness of each pixel:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "digits = datasets.load_digits()\n",
    "digits.data.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's show a quick classification example, using the simple-and-fast Gaussian Naive Bayes estimator."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print (X)\n",
    "print (X[:-5])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.naive_bayes import GaussianNB\n",
    "X = digits.data\n",
    "y = digits.target\n",
    "\n",
    "# Instantiate the estimator\n",
    "clf = GaussianNB()\n",
    "\n",
    "# Fit the estimator to the data, leaving out the last five samples\n",
    "clf.fit(X[:-5], y[:-5])\n",
    "\n",
    "# Use the model to predict the last several labels\n",
    "y_pred = clf.predict(X[-5:])\n",
    "\n",
    "print (y_pred)\n",
    "print (y[-5:])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We see that this relatively simple model leads to a perfect classification of the last few digits!\n",
    "\n",
    "Let's use the model to predict labes for the full dataset, and plot the **confusion matrix**, which is a convenient visual representation of how well the classifier performs:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn import metrics\n",
    "\n",
    "clf = GaussianNB()\n",
    "clf.fit(X, y)\n",
    "y_pred = clf.predict(X)\n",
    "\n",
    "def plot_confusion_matrix(y_pred, y):\n",
    "    plt.imshow(metrics.confusion_matrix(y, y_pred),\n",
    "               cmap=plt.cm.binary, interpolation='nearest')\n",
    "    plt.colorbar()\n",
    "    plt.xlabel('true value')\n",
    "    plt.ylabel('predicted value')\n",
    "    \n",
    "print (\"classification accuracy:\", metrics.accuracy_score(y, y_pred))\n",
    "plot_confusion_matrix(y, y_pred)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Interestingly, there is confusion between some values.  In particular, the number **2** is often mistaken for the number **8** by this model!  But for the vast majority of digits, we can see that the lassification looks correct.\n",
    "\n",
    "Let's use the ``metrics`` submodule again to print the accuracy of the classification:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print (metrics.accuracy_score(y, y_pred))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We have an 82% accuracy rate with this particular model.\n",
    "\n",
    "But there's a problem: we are testing the model on the data we used to train the model. As we'll see later, this is generally not a good approach to model validation!  Because of the nature of the Naive Bayes estimator, it's alright in this case, but we'll see later examples where this approach causes problems."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Measuring Prediction Performance\n",
    "\n",
    "An important piece of the learning task is the measurement of prediction performance, also known as **model validation**.  We'll go into detail about this, but first motivate the approach with an example.\n",
    "\n",
    "## The Importance of Splitting\n",
    "Above we looked at a *confusion matrix*, which can be computed based on the results of any model. Let's look at another classification scheme here, the *K-Neighbors Classifier*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.neighbors import KNeighborsClassifier\n",
    "from sklearn import datasets\n",
    "\n",
    "digits = datasets.load_digits()\n",
    "X, y = digits.data, digits.target\n",
    "\n",
    "clf = KNeighborsClassifier(n_neighbors=1)\n",
    "clf.fit(X, y)\n",
    "y_pred = clf.predict(X)\n",
    "\n",
    "print (\"classification accuracy:\", metrics.accuracy_score(y, y_pred))\n",
    "plot_confusion_matrix(y, y_pred)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Our classifier gives perfect results!  Have we settled on a perfect classification scheme?\n",
    "\n",
    "**No!**  The *K*-neighbors classifier is an example of an instance-based classifier, which memorizes the input data and compares any unknown sample to it.  To accurately measure the performance, we need to use a separate *validation set*, which the model has not yet seen.\n",
    "\n",
    "Scikit-learn contains utilities to split data into a training and validation set:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.cross_validation import train_test_split\n",
    "\n",
    "X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3)\n",
    "print (X_train.shape, X_test.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "clf = KNeighborsClassifier(n_neighbors=1)\n",
    "clf.fit(X_train, y_train)\n",
    "y_pred = clf.predict(X_test)\n",
    "\n",
    "print (\"classification accuracy:\", metrics.accuracy_score(y_test, y_pred))\n",
    "plot_confusion_matrix(y_test, y_pred)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This gives us a more accurate indication of how well the model is performing.\n",
    "\n",
    "For this reason you should **always do a train/test split** when validating a model."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exploring Validation Metrics\n",
    "\n",
    "Above, we used perhaps the most simple evaluation metric, the number of matches and mis-matches.  But this is not always sufficient.  For example, imagine you have a situation where you'd like to identify a rare class of event from within a large number of background sources."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Generate an un-balanced 2D dataset\n",
    "np.random.seed(0)\n",
    "X = np.vstack([np.random.normal(0, 1, (950, 2)),\n",
    "               np.random.normal(-1.8, 0.8, (50, 2))])\n",
    "y = np.hstack([np.zeros(950), np.ones(50)])\n",
    "\n",
    "plt.scatter(X[:, 0], X[:, 1], c=y, edgecolors='none',\n",
    "            cmap=plt.cm.Accent)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exploring other Validation Scores\n",
    "\n",
    "Until now we are using only the **accuracy** to evaluate our algorithms. We can calculate other scores such as the **precision**, the **recall**, and the **f1 score**:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn import metrics\n",
    "from sklearn.svm import SVC\n",
    "\n",
    "X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=0)\n",
    "clf = SVC(kernel='linear').fit(X_train, y_train)\n",
    "y_pred = clf.predict(X_test)\n",
    "\n",
    "print (\"accuracy:\", metrics.accuracy_score(y_test, y_pred))\n",
    "print (\"precision:\", metrics.precision_score(y_test, y_pred))\n",
    "print (\"recall:\", metrics.recall_score(y_test, y_pred))\n",
    "print (\"f1 score:\", metrics.f1_score(y_test, y_pred))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### What do these mean?\n",
    "\n",
    "These are ways of taking into account not just the classification results, but the results **relative to the true category**.\n",
    "\n",
    "- $$ {\\rm accuracy} \\equiv \\frac{\\rm correct~labels}{\\rm total~samples} $$\n",
    "\n",
    "- $$ {\\rm precision} \\equiv \\frac{\\rm true~positives}{\\rm true~positives + false~positives} $$\n",
    "\n",
    "- $$ {\\rm recall} \\equiv \\frac{\\rm true~positives}{\\rm true~positives + false~negatives} $$\n",
    "\n",
    "- $$ F_1 \\equiv 2 \\frac{\\rm precision \\cdot recall}{\\rm precision + recall} $$\n",
    "\n",
    "The **accuracy**, **precision**, **recall**, and **f1-score** all range from 0 to 1, with 1 being optimal.\n",
    "Here we've used the following definitions:\n",
    "\n",
    "- *True Positives* are those which are labeled ``1`` which are actually ``1``\n",
    "- *False Positives* are those which are labeled ``1`` which are actually ``0``\n",
    "- *True Negatives* are those which are labeled ``0`` which are actually ``0``\n",
    "- *False Negatives* are those which are labeled ``0`` which are actually ``1``\n",
    "\n",
    "\n",
    "We can quickly compute a summary of these statistics using scikit-learn's provided convenience function:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print (metrics.classification_report(y_test, y_pred,\n",
    "                                    target_names=['background', 'foreground']))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This tells us that, though the overall correct classification rate is 97%, we only correctly identify 67% of the desired samples, and those that we label as positives are only 83% correct!  This is why you should make sure to carefully choose your metric when validating a model."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Cross-Validation\n",
    "\n",
    "Using the simple train/test split as above can be useful, but there is a disadvantage: **You're ignoring a portion of your dataset**.  One way to address this is to use cross-validation.\n",
    "\n",
    "The simplest cross-validation scheme involves running two trials, where you split the data into two parts, first training on one, then training on the other:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X1, X2, y1, y2 = train_test_split(X, y, test_size=0.5)\n",
    "print (X1.shape)\n",
    "print (X2.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "y2_pred = SVC(kernel='linear').fit(X1, y1).predict(X2)\n",
    "y1_pred = SVC(kernel='linear').fit(X2, y2).predict(X1)\n",
    "\n",
    "print (np.mean([metrics.precision_score(y1, y1_pred),\n",
    "               metrics.precision_score(y2, y2_pred)]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This is known as **two-fold** cross-validation, and is a special case of *K*-fold cross validation.\n",
    "\n",
    "Because it's such a common routine, scikit-learn has a K-fold cross-validation scheme built-in:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.cross_validation import cross_val_score\n",
    "\n",
    "# Let's do a 2-fold cross-validation of the SVC estimator\n",
    "print (cross_val_score(SVC(kernel='linear'), X, y, cv=2, scoring='precision'))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It's also possible to use ``sklearn.cross_validation.KFold`` and ``sklearn.cross_validation.StratifiedKFold`` directly, as well as other cross-validation models which you can find in the ``cross_validation`` module."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Example: The ``SVC`` classifier takes a parameter ``C`` whose default value is ``1``.  Using 5-fold cross-validation, make a plot of the precision as a function of ``C``, for the ``SVC`` estimator on this dataset.  For best results, use a logarithmic spacing of ``C`` between 0.01 and 100."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Cs = logspace(-2, 2, 40)\n",
    "scores = []\n",
    "\n",
    "for C in Cs:\n",
    "    score = cross_val_score(SVC(kernel='linear', C=C), X, y, cv=5, scoring='precision')\n",
    "    scores.append(score.mean())\n",
    "\n",
    "plot(Cs, scores)\n",
    "xscale('log')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Grid Search\n",
    "\n",
    "The previous exercise is an example of a **grid search** for model evaluation.  Again, because this is such a common task, Scikit-learn has a grid search tool built-in, which is used as follows.  Note that ``GridSearchCV`` has a ``fit`` method: it is a meta-estimator: an estimator over estimators!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.grid_search import GridSearchCV\n",
    "\n",
    "clf = SVC(kernel='linear')\n",
    "Crange = np.logspace(-2, 2, 40)\n",
    "\n",
    "grid = GridSearchCV(clf, param_grid={'C': Crange},\n",
    "                    scoring='precision', cv=5)\n",
    "grid.fit(X, y)\n",
    "\n",
    "print (\"best parameter choice:\", grid.best_params_)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "scores = [g[1] for g in grid.grid_scores_]\n",
    "plt.semilogx(Crange, scores);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Grid search can come in very handy when you're tuning a model for a particular task."
   ]
  }
 ],
 "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": 1
}