5 - Modeling with scikit-learn.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Modeling with Scikit-Learn"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this notebook, we're going to tie everything we've learned so far to do some modeling with scikit-learn. As we've seen, we have in the violation descriptions a large number of, relatively, free text fields.\n",
    "\n",
    "Working with text data is a particularly attractive use case for machine learning. It's also often a messy one that can involve working with a lot of boilerplate code. The library [scikit-learn](http://scikit-learn.org/stable/) provides many features for working with text data.\n",
    "\n",
    "First, let's take a closer look at scikit-learn."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Preliminaries: Scikit-Learn\n",
    "\n",
    "The `scikit-learn` package provides a robust set of machine learning algorithms for Python. Like all of the packages, we have seen so far, scikit-learn is built upon the core Python scientific stack (i.e. NumPy, SciPy, Cython). One of the biggest reasons of why scikit-learn is so popular is that it has a simple, consistent API, making it useful for a wide range of statistical learning applications. The different components of scikit-learn can be combined to make powerful and expressive pipelines for analyzing data.\n",
    "\n",
    "Scikit-learn provides facilities for\n",
    "\n",
    "* **supervised learning** algorithms that learn from a training set with **labels**, or targets, to generalize to other inputs like **regression** and **classification**.\n",
    "* **unsupervised learning** algorithms that learn structure in the data from a training set of unlabeled examples like **clustering** or **density estimators**\n",
    "* **dimensionality reduction** algorithms which reduce the number of **features**, or columns, while preserving information about the data\n",
    "* **model selection** for choosing the best parameters and models\n",
    "* **preprocessing** for getting data ready to apply machine learning algorithms"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Representing Data in Scikit-Learn\n",
    "\n",
    "Most machine learning algorithms implemented in scikit-learn expect data to be stored in a two-dimensional array or matrix. The arrays can be either **numpy** arrays, or in some cases **scipy.sparse** matrices. 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). A sample can be a document, a picture, a sound, a video, a row in database or CSV file, or whatever you can describe with a fixed set of quantitative traits.\n",
    "\n",
    "* **n_features**: The number of features or distinct traits that can be used to describe each item in a quantitative manner. Features are generally real-valued, but may be boolean or discrete-valued in some cases.\n",
    "\n",
    "The number of features (almost always) must be fixed in advance. However it can be very high dimensional (e.g. millions of features) with most of them being zeros for a given sample. This is a case where `scipy.sparse` matrices and other techniques can be useful, in that they are much more memory-efficient than numpy arrays."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Aside: NumPy Arrays\n",
    "\n",
    "We haven't talked much about NumPy arrays. NumPy arrays, however, are the fundamental data structure in Python data stack. \n",
    "\n",
    "A NumPy array is an object that represents a homogeneously typed, multidimensional array. The array provides an efficient (close to the hardware) data structure for scientific, or array-oriented, computing. First, let's look at the NumPy import convention."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "import numpy as np"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And create an array from a Python list. This is an array of all integers."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "x = np.array([1, 2, 3, 4, 5])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can't assign something that's not an integer."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "x[0] = 'a'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can perform indexing operations much like we saw with pandas earlier, but without the convenience of labels."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can use regular Python slicing syntax."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "x[:3]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "x[::2]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Or what's called **fancy indexing** by using Boolean or integer indexes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "x[[True, False, True, False, True]]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can perform operations on NumPy arrays like `sum`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "np.array([1, 2, 3, 4, 5]).sum()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And we can perform linear algebra operations, like taking the dot product."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "x = np.array([[1, 2, 3], \n",
    "              [4, 5, 6],\n",
    "              [7, 8, 9]])\n",
    "\n",
    "x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "y = np.array([[4], [5], [6]])\n",
    "\n",
    "y"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "x.dot(y)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Array Creation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Both the `zeros` and `ones` functions can be useful for creating data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.zeros(5, dtype=float)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.zeros(5, dtype=int)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.zeros(5, dtype=complex)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.ones(5, dtype=float)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Two other handy functions to know about are `arange`, `linspace`, and `logspace`.\n",
    "\n",
    "`np.arange` creates an array of a range of integers."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.arange(10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`np.linspace` and `np.logspace` create linearly and logarithmically-spaced grids, respectively, with a fixed number of points."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.linspace(0, 1, num=5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.logspace(-1, 3, num=5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise\n",
    "\n",
    "Create an array with 1000 numbers, 0 to 999, using `arange`. Have a look at the `reshape` method on arrays to turn this into an array with 10 columns."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Type your solution here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%load solutions/numpy_reshape.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, it is often useful to create arrays with random numbers that follow a specific distribution. The np.random module contains a number of functions that can be used to this effect, for example this will produce an array of 5 random samples taken from a standard normal distribution (0 mean and variance 1) $ X \\sim N(0, 1) $:\n",
    "\n",
    "$$f(x \\mid \\mu = 0, \\sigma=1) = \\sqrt{\\frac{1}{2\\pi\\sigma^2}}\\exp {-\\frac{x^2}{2\\sigma^2} }$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.random.randn(5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$X \\sim N(9, 3)$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "norm10 = np.random.normal(loc=9, scale=3, size=10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Generate a NumPy array of 1000 random numbers sampled from a Poisson distribution, with parameter `lam=5`. What is the modal value in the sample? You maybe interested in using `np.bincounts`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Type your solution here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%load solutions/random_poisson.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "NumPy and the SciPy libraries also provide much more than data structures like more facilities for linear algebra, matrix decompositions, optimization, clustering, polynomials, unit testing, etc."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Scikit-Learn Quickstart"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's take a quick look at scikit-learn to fix ideas before going much further. We'll have a look at the canonical iris dataset, which consists of a set of measurements for flowers, each being a member of one of three species: Iris Setosa, Iris Versicolor or Iris Virginica."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "from sklearn.datasets import load_iris"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "dta = load_iris()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The features of the data consists of"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "dta.feature_names"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The labels consist of"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "dta.target_names"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "dta.data[:10]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "dta.target[::10]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's fit a logistic regression model on all of the iris data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "from sklearn.linear_model import LogisticRegression"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "model = LogisticRegression()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "model.fit(dta.data, dta.target)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "model.predict(dta.data)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### `scikit-learn` interface\n",
    "\n",
    "The power of scikit-learn comes from the fact that they share a common, unified API, consisting of three complementary interfaces:\n",
    "\n",
    "* **estimator** interface for building and ﬁtting models\n",
    "* **predictor** interface for making predictions\n",
    "* **transformer** interface for converting data.\n",
    "\n",
    "The estimator interface is at the core of the library. It deﬁnes instantiation mechanisms of objects and exposes a fit method for learning a model from training data. All supervised and unsupervised learning algorithms (*e.g.*, for classiﬁcation, regression or clustering) are oﬀered as objects implementing this interface. Machine learning tasks like feature extraction, feature selection or dimensionality reduction are also provided as estimators.\n",
    "Scikit-learn strives to have a uniform interface across all methods. For example, a typical **estimator** follows this template:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```python\n",
    "class Estimator:\n",
    "  \n",
    "    def fit(self, X, y=None):\n",
    "        \"\"\"Fit model to data X (and y)\"\"\"\n",
    "        self.some_attribute = self.some_fitting_method(X, y)\n",
    "        return self\n",
    "            \n",
    "    def predict(self, X_test):\n",
    "        \"\"\"Make prediction based on passed features\"\"\"\n",
    "        pred = self.make_prediction(X_test)\n",
    "        return pred\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For a given scikit-learn estimator object named model, several methods are available. Irrespective of the type of estimator, there will be a fit method:\n",
    "\n",
    "* model.fit : fit training data. For supervised learning applications, this accepts two arguments: the data X and the labels y (e.g. `model.fit(X, y)`). For unsupervised learning applications, this accepts only a single argument, the data X (e.g. `model.fit(X)`).\n",
    "\n",
    "> During the fitting process, the state of the **estimator** is stored in attributes of the estimator instance named with a trailing underscore character (`_`). For example, the sequence of regression trees `sklearn.tree.DecisionTreeRegressor` is stored in `estimators_` attribute."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The **predictor** interface extends the notion of an estimator by adding a predict method that takes an array `X_test` and produces predictions based on the learned parameters of the estimator. In the case of supervised learning estimators, this method typically returns the predicted labels or values computed by the model. Some unsupervised learning estimators may also implement the predict interface, such as k-means, where the predicted values are the cluster labels."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "all **supervised estimators** are expected to have the following methods:\n",
    "\n",
    "* `model.predict` : given a trained model, predict the label of a new set of data. This method accepts one argument, the new data X_new (e.g. `model.predict(X_new)`), and returns the learned label for each object in the array.\n",
    "* `model.predict_proba` : For classification problems, some estimators also provide this method, which returns the probability that a new observation has each categorical label. 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 a score method. Scores are between 0 and 1, with a larger score indicating a better fit."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since it is common to modify or ﬁlter data before feeding it to a learning algorithm, some estimators in the library implement a **transformer** interface which deﬁnes a `transform` method. It takes as input some new data `X_test` and yields as output a transformed version. Preprocessing, feature selection, feature extraction and dimensionality reduction algorithms are all provided as transformers within the library."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**unsupervised estimators** will always have these methods:\n",
    "\n",
    "* `model.transform` : given an unsupervised model, transform new data into the new basis. This also accepts one argument  `X_new`, and returns the new representation of the data based on the unsupervised model.\n",
    "* `model.fit_transform` : some estimators implement this method, which more efficiently performs a fit and a transform on the same input data."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's take a look at some examples of each of these using the Chicago Health Inspection data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "from load_data import dta"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dta = dta.drop([\"violations\"], axis='columns').join(\n",
    "    dta.violations.str.split(\"|\", expand=True)\n",
    "    .unstack()\n",
    "    .dropna()\n",
    "    .str.strip()\n",
    "    .reset_index(level=0, drop=True)\n",
    "    .to_frame()\n",
    "    .rename(columns={0: 'violations'}),\n",
    "    how='right'\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dta.info()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Bag-of-Words\n",
    "\n",
    "First, we need to take our text and turn it in to numerical features. A common assumption for doing machine learning on text is what's known as the bag of words assumption. This means that we assume that the order of the words as they occur in a document doesn't matter to discern the general meaning of the document. This is commonly done in the following steps\n",
    "\n",
    "1. Build what's called a **vocabulary**, which is a mapping from integers to possible words, $w$, in your corpus, or collection of documents.\n",
    "2. Using this vocabulary, assign a number to the count of each word occuring in any document.\n",
    "\n",
    "What you're left with is a matrix $X$, where each value $X[i,j]$ is the count of word $j$ in document $i$.\n",
    "$X$ is a matrix of dimension `n_documents` by `n_vocabulary`. This is large. Luckily, most words don't occur in every document. If they did, we would not be able to separate the documents according to topics.\n",
    "\n",
    "For this reason, bag of words documents are often high-dimensional, sparse datasets. We don't need to keep the zeros in memory."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Tokenize"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Ok, so how do we do this? Text is often really messy, has punctuation, and has a bunch of words that every text has to have but don't necessarily connote topical meaning. These words are called stop words such as \"the,\" \"a,\" or \"an.\"\n",
    "\n",
    "We turn human writing into a set of feature vectors by taking care of these issues. This process is called tokenization.\n",
    "\n",
    "scikit-learn provides some nice facilities for building a dictionary of features and transforming documents to feature vectors. The first of these that we will look at is the **CountVectorizer** transformer.\n",
    "\n",
    "Recall from above that a transformer is an estimator that provides a transform method."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "from sklearn.feature_extraction.text import CountVectorizer"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In scikit-learn, all of the estimators take their options when you instantiate the estimator. Here, we say that we want to remove stop-words using a list of common english language words that we won't need."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "count_vectorizer = CountVectorizer(stop_words='english')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then we need to *fit* the transformer on the data. Calling fit always returns the object itself. We'll see why later.\n",
    "\n",
    "Calling a `fit` method on an estimator actually does the learning. Any learned parameters are now attached to the estimator with an underscore (`_`) appended."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "count_vectorizer.fit(dta.violations)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the case of `CountVectorizer` there is a dictionary called `vocabulary_` which stores a mapping from the known vocabulary to the column in the sparse matrix which contains the counts for that word. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "len(count_vectorizer.vocabulary_)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, we need to transform our original data using transform."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "count_matrix = count_vectorizer.transform(dta.violations)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Count matrix is a **sparse matrix**, provided by the SciPy library. The number of samples is equal to the number of violations that we have in the data. The number of columns is the cardinality of our vocabulary. The entries are the counts of each word in the document."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "count_matrix"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Sparse matrices behave a lot like plain numpy arrays. For example, we can ask for the sum of each word over all the documents."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "count_matrix.sum(0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise\n",
    "\n",
    "What is the most frequent word in this vocabulary? Explore the `count_vectorizer` method to see if it offers anything helpful in uncovering this."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Type your solution here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "%load solutions/frequent_word.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Tf-Idf"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Looking at the most common word we already see one issue with using raw counts. Another issue is that longer documents will have higher counts of words. Commonly, we use a technique known as **term frequency - inverse document frequency**, or **tf-idf**, instead of counts to do analysis on text data, which mitigates these issues.\n",
    "\n",
    "The *term frequency* is a measure of the frequency of a word in a document. Term frequency in document $i$ for word $j$ is\n",
    "\n",
    "$$tf_{ij}=\\frac{w_{ij}}{\\sum_jw_{ij}}$$\n",
    "\n",
    "You might go about computing this."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.preprocessing import normalize\n",
    "\n",
    "tf = normalize(count_matrix, norm='l1', axis=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using the `l1`-norm across axis 1 (the column) index, we now have frequencies within the document. Each document now sums to 1.\n",
    "\n",
    "One thing to point out here is that summations over a `scipy.sparse` matrix returns a numpy `matrix`. This is mostly for historical reasons, and I *don't* recommend working with the `matrix` data structure if you can avoid it. You can turn this into an array by accessing the matrix's `A` attribute."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "tf.sum(1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "tf.sum(1).A"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Another important concept is that of inverse document frequency. This is a measure of how important a word is. Words like stop words or words that are otherwise popular in a corpus will still have a high term frequency. Inverse document frequency is a way to downweight the frequent terms but upweight the rare ones. The inverse document frequency is\n",
    "\n",
    "$$idf = \\log\\left(\\frac{N_{\\text{documents}}}{N_{\\text{documents with term}}}\\right)$$\n",
    "\n",
    "or\n",
    "\n",
    "$$idf = \\log\\left(\\frac{N_{\\text{documents}}}{1 + N_{\\text{documents with term}}}\\right)$$\n",
    "\n",
    "in case your vocabulary is a superset of the words in your documents.\n",
    "\n",
    "So tf-idf is\n",
    "\n",
    "$$\\text{tf-idf} = tf \\times idf$$\n",
    "\n",
    "Scikit-learn actually uses a slightly different definition.\n",
    "\n",
    "Of course, scikit-learn provides a transformer for tf-idf"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "from sklearn.feature_extraction.text import TfidfVectorizer"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's prepare our TfidfVectorizer. We'll remove stop-words, remove any words that don't occur in at least 50 documents and remove words that occur in 85% or more documents.\n",
    "\n",
    "Finally, we'll use a **regular expression** pattern to determine what exactly a token (or word) is. In this case, we deviate from the scikit-learn default by not allowing numbers to be words."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "tfidf_vect = TfidfVectorizer(\n",
    "    stop_words='english', \n",
    "    min_df=50,\n",
    "    max_df=.85, \n",
    "    token_pattern=r\"(?u)\\b[A-Za-z_][A-Za-z_]+\\b\"\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Notice here that we can combine the fitting and the transformation by taking advantage of the `fit_transform` method."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "X = tfidf_vect.fit_transform(dta.violations)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "X"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using these more restrictive criteria above, we've greatly reduced the dimensionality of the feature space, while ideally preserving the most useful information on the contents of the documents."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Dimensionality Reduction"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's take a look at another kind of transformer in scikit-learn, one that provides dimensionality reduction. Here we'll use Truncated SVD on the tf-idf matrix. Formally, this is known as Latent Semantic Analysis (LSA), because it transforms the documents to a low-dimensional \"semantic\" space. Truncated SVD is a lot like Principle Components Analysis (PCA), except that the decomposition is on the documents rather than the covariance matrix. \n",
    "\n",
    "\n",
    "Mathematically, truncated SVD applied to training samples X produces a low-rank approximation $X_k$:\n",
    "\n",
    "$$X \\approx X = U_k \\Sigma_k V_k^\\top$$\n",
    "\n",
    "After this operation, $U_k \\Sigma_k^\\top$ is the transformed training set with k features (called `n_components` in the API).\n",
    "\n",
    "To also transform a test set $X$, we multiply it with $V_k$:\n",
    "\n",
    "$$X' = X V_k$$\n",
    "\n",
    "If we were to center the matrix $X$ then TruncatedSVD would be equivalent to PCA. Not doing so allows us to continue to work with sparse matrices as documents almost always produce."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`TruncatedSVD` is available under the `decomposition` namespace."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "from sklearn.decomposition import TruncatedSVD"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here we'll fit the `TruncatedSVD` transformer using 10 components. This is an arbitrary choice. In practice, you may want to tune the number of components using whatever metric is appropriate for your task. Note that we use `random_state` here to make sure that our results are repeatable. Any of the algorithms in scikit-learn that are non-deterministic will provide a `random_state` keyword. It is really important that you use it to ensure **repeatable** results."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "n_components = 10\n",
    "\n",
    "svd = TruncatedSVD(\n",
    "    n_components=n_components, \n",
    "    random_state=0\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We'll use `fit_transform` to perform the singular value decomposition (up to the first $k$ components) and to project the original matrix into the reduced space."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "X_reduced = svd.fit_transform(X)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_reduced"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_reduced.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "words = np.array(sorted(tfidf_vect.vocabulary_.keys()))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "words[:15]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise\n",
    "\n",
    "Write a loop that prints the top ~6 words for each component according to the the magnitude of their loadings."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Type your solution here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "%load solutions/top_words_loadings.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Clustering"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now that we've done some dimensionality reduction, we may be interested in clustering the documents in this reduced space. Scikit-Learn has a number of clustering algorithms. Here we'll use K-means.\n",
    "\n",
    "The K-means algorithm clusters data by separating it into groups of equal variance, choosing them in order to minimize the within-cluster sum of squares. Formally, we divide $n$ samples into $k$ clusters $C$. Each cluster is defined by its mean, or centroid, $u_i$.\n",
    "\n",
    "$$\\sum_{i=0}^n\\underset{u_j\\in C}\\min (\\|x_j - u_i\\|^2)$$\n",
    "\n",
    "K-means proceeds as follows:\n",
    "\n",
    "1. We pick $k$ random points from the dataset and call them the cluster centroids\n",
    "2. We assign each data point to its closest centroid.\n",
    "3. We recompute the centroids.\n",
    "4. The distance between the old and new centroids are computed until they stop moving.\n",
    "\n",
    "Eventually k-means will converge. However, there is no guarantee that it will converge to a global optimum. One thing we can do to mitigate this is to pick better starting points than $k$ random points in the data. Scikit-learn uses a better choice by default through the `init='k-means++'` argument, which attempts to pick starting centroids that are generally 'far' from each other.\n",
    "\n",
    "First, let's normalize the data row-wise. If we do this, for documents the euclidean distance above becomes cosine similarity. Now we're performing spherical k-means so that all of our comparisons between documents are equal and indpendent of the size of the document."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "from sklearn.preprocessing import Normalizer\n",
    "\n",
    "normalizer = Normalizer(copy=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "X_norm = normalizer.fit_transform(X_reduced)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "np.linalg.norm(X_norm, axis=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here, we have some guidance on the number of clusters, we might expect, so let's use it.\n",
    "\n",
    "Recall that we found 45 distinct violation numbers in the data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n_clusters = (dta.violations.str.extract(\"(\\d+)(?=\\.)\", expand=False)\n",
    "              .astype(int).unique().shape[0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n_clusters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "from sklearn.cluster import KMeans\n",
    "\n",
    "kmeans = KMeans(n_clusters=n_clusters, random_state=0)\n",
    "\n",
    "kmeans.fit(X_norm)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise\n",
    "\n",
    "Plot the histogram of the number of found clusters vs. our known violations. They won't line up exactly. I.e., the bins found for clustering will be different than the extracted violation numbers.\n",
    "\n",
    "What are the first five violations in, say, the first three clusters?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "# Type your solution here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%load solution/plot_clusters.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise\n",
    "\n",
    "See if you can strip out the comments and still find (semi-)meaningful clusters. Here's a hint, you'll want to again use regular expressions and the `str` accessor in pandas."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here, we can use a *lookbehind* (`(?<=)`) to capture via `()` one or more (`+`) of any character (`.`) that follows the word \"Comments.\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import re\n",
    "\n",
    "result = re.search(\"(?<=Comments:)(.+)\", \n",
    "                   \"1. This is a violation. Comments: This was a really egregious violation.\")\n",
    "\n",
    "result"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "result.group()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Type your solution here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%load solutions/clustering_comments.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Modeling"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Can we predict a pass/fail rating from features in the data? First, let's start to build our modeling set by turning the data back into inspection level data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dta.columns"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For the sake of clarity, let's drop all inspection results that weren't Pass or Fail. We can use another DataFrame method `isin` to do this."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dta.results.unique()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dta = dta.loc[dta.results.isin(['Pass', 'Fail'])]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let's pull out the columns that we think may be helpful in predicting."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "columns = [\n",
    "    \"inspection_id\",\n",
    "    \"inspection_date\", \n",
    "    \"inspection_type\", \n",
    "    \"facility_type\", \n",
    "    \"results\", \n",
    "    \"risk\",\n",
    "    \"zip\"\n",
    "]\n",
    "\n",
    "modeling_dta = (\n",
    "    dta.reset_index()\n",
    "    .drop_duplicates([\"inspection_id\"])\n",
    "    .loc[:, columns]\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "modeling_dta.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's do some **feature engineering** to turn our data into some more features we think may be predictive of passing or failing scores."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We might use the number of violations per establishment during an inspection."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "number_of_violations = dta.groupby(dta.index).size()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Maybe some violation severity is weather dependent (on average)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "modeling_dta[\"month\"] = modeling_dta.inspection_date.dt.month.astype(str)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, let's add in the number of previously failed inspections. We'll ignore the fact that we may care more about the rate of failed inspections or the full history of failures."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "key = (\n",
    "    \"address\", \n",
    "    \"dba_name\", \n",
    ")\n",
    "\n",
    "\n",
    "def failed_last(df):\n",
    "    # return a Series to name the result and preserve the index\n",
    "    # it's important to use `values` here or pandas will try to\n",
    "    # align the indices you provide vs. what's in the series\n",
    "    return pd.Series((df.shift(1).results == 'Fail').cumsum().values, \n",
    "                     name='num_fails',\n",
    "                     index=df.inspection_id)\n",
    "\n",
    "\n",
    "fail_num = (dta\n",
    "            .reset_index()\n",
    "            .drop_duplicates(['inspection_id'])\n",
    "            .sort_values(['address', 'dba_name', 'inspection_date'])\n",
    "            .groupby(key)\n",
    "            .apply(failed_last)\n",
    "            .reset_index(level=[0, 1], drop=True)\n",
    "            .reset_index())\n",
    "\n",
    "fail_num.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "modeling_dta = modeling_dta.merge(fail_num)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "I always like to check shapes after merging to make sure I got the results I expected. By the way, there is a fantastic library [engarde](https://github.com/TomAugspurger/engarde) for doing direct defensive/assertive programming. It allows you to write functions like\n",
    "\n",
    "\n",
    "```python\n",
    "@is_shape((1290, 10))\n",
    "@unique_index\n",
    "def make_design_matrix('data.csv'):\n",
    "    out = ...\n",
    "    return out\n",
    "```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "modeling_dta.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "modeling_dta.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's join the violations data back together and run truncated SVD on the tf-idf matrix to reduce the dimensionality."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "violations = dta.violations.groupby(dta.index).apply(lambda df: \" \".join(df))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Again, let's be defensive here and make sure our violations data is sorted the same as our modeling data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "violations = violations.loc[modeling_dta.set_index('inspection_id').index]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X = tfidf_vect.fit_transform(violations)\n",
    "\n",
    "svd = TruncatedSVD(n_components=10, random_state=0)\n",
    "\n",
    "X_reduced = svd.fit_transform(X)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The last thing we'll need to do is to turn all of these non-numeric features into numeric features."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "modeling_dta.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dummied = pd.get_dummies(\n",
    "    modeling_dta[['inspection_type', 'facility_type', 'risk', 'zip', 'month', 'num_fails']]\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dummied.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*Technically* you probably don't want to do it this way. We'll discuss why below."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can now join these values with the LSA results above."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X = np.column_stack((dummied, X_reduced))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's pull out the field that we want to try to predict, `results`. We'll use another scikit-learn transformer, `LabelEncoder` to do prepare the Pass/Fail column for modeling.\n",
    "\n",
    "The `LabelEncoder` learns and transforms any labels to values from `0` to `n_classes - 1`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.preprocessing import LabelEncoder"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "encoder = LabelEncoder()\n",
    "\n",
    "y = encoder.fit_transform(modeling_dta.results)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "y"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since the `LabelEncoder` preserves lexicographical ordering, Pass is now 1."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "When we're doing predictive modeling, the main thing we're concerned with is 'does my model generalize.' Will it accurately predict on data that the model hasn't seen during training.\n",
    "\n",
    "To answer this question, we might split our dataset into a sample to train on and a separate sample to test on. Let's again use scikit-learn to prepare our train and test sets."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.model_selection import train_test_split\n",
    "\n",
    "(X_train, X_test, \n",
    " y_train, y_test) = train_test_split(X, y, test_size=.25, random_state=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's check that our classes were relatively evenly distributed across our train and test splits. There are better ways to ensure this with scikit-learn, but we will just do a quick spot check."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "y_train.mean()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "y_test.mean()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Aside: Decision Trees\n",
    "\n",
    "A decision tree is a model that recursively partitions a dataset in such a way as to lead to the largest information gain. For each partition, a simple model is fit, say, a constant. Decision trees are good at capturing non-linearities and feature interactions by design.\n",
    "\n",
    "Here we'll look at the `DecisionTreeClassifier` model estimator from scikit-learn."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.tree import DecisionTreeClassifier, export_graphviz\n",
    "\n",
    "dtree = DecisionTreeClassifier(random_state=0, max_depth=4)\n",
    "\n",
    "dtree.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "columns = dummied.columns.tolist()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n_components = X_reduced.shape[1]\n",
    "\n",
    "columns += ['x{}'.format(i) for i in range(n_components)]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "columns"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Note**: If you can't run this part, don't worry. You'll need to have graphviz installed. (`conda install graphviz` *should* do the trick.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "export_graphviz(dtree, out_file='tree.dot', feature_names=columns)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!dot -Tpng tree.dot -o tree.png"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from IPython.display import Image\n",
    "Image('tree.png')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "How well does our decision tree do out of sample? We can use `predict_proba` method available on most classifiers to get our predicted probabilities for each sample for each class."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "proba = dtree.predict_proba(X_test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "proba"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We might check the area under the receiver operating characteristic curve to see how well our shallow tree does. Intuitively, this measure tells you if you were to grab a sample at random, how likely is your classifier to predict the correct label. Scikit-learn provides a number of useful metrics under the `sklearn.metrics` namespace."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.metrics import roc_auc_score\n",
    "\n",
    "roc_auc_score(y_test, proba[:, 1])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we'll see if we can do slightly better using a RandomForest a model. The RandomForest model is, roughly, an ensemble of decision trees. We will fit 500 full-depth trees, selecting features and rows at random, and then average the votes from each tree to get our classification results."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.ensemble import RandomForestClassifier"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "clf = RandomForestClassifier(n_estimators=500, n_jobs=4)\n",
    "\n",
    "clf"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "clf.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "All of the model estimators have a `score` method, which is by default the average accuracy."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "clf.score(X_test, y_test)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Of course, we have a roughly 1:3 class imbalance here, so we may also want to check the AUC score and have a look at the confusion matrix.\n",
    "\n",
    "In the confusion matrix the true labels are the rows, and the predicted labels are the columns. Here we see that we're slightly high on our false positive rate. Inspections that were a Fail are being predicted as Pass."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.metrics import confusion_matrix"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "confusion_matrix(y_test, clf.predict(X_test))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "roc_auc_score(y_test, clf.predict_proba(X_test)[:, 1])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's make `columns` into a numpy array, so we can use fancy-indexing."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "columns = np.array(columns)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using the `feature_importances_` attribute of the trained Random Forest, we can get a sense of what features are important. Roughly speaking this is a measure across the decision trees of how influential a feature was in improving predictions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "clf.feature_importances_"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's rank the features in descending order and see which ones are the most important."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "idx = clf.feature_importances_.argsort()[::-1]\n",
    "\n",
    "list(zip(clf.feature_importances_[idx][:35], columns[idx][:35]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Putting It All Together"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Consider all of the steps that we've been through so far. We've\n",
    "\n",
    "* loaded our data\n",
    "* cleaned our data\n",
    "* extracted and engineered features\n",
    "* fit a model\n",
    "* evaluated this model, and\n",
    "* iterated on these steps\n",
    "\n",
    "Scikit-learn is designed from the ground up to make these steps easy for users with minimal boilerplate. Recall the main scikit-learn objects -- the `Transformer`, `Predictor`, and `Estimator`.\n",
    "\n",
    "The scikit-learn `Pipeline` abstraction builds on these interfaces to allow us to put together a chain of transformers and estimators and use the pipeline, as if it were an estimator itself."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Transformer Interface\n",
    "\n",
    "Recall the transformer interface. A transformer is intended to filter or modify the data in a supervised or unsupervised way.\n",
    "\n",
    "```python\n",
    "new_data = obj.transform(data)\n",
    "```\n",
    "\n",
    "The interface is\n",
    "\n",
    "```python\n",
    "class Transformer:\n",
    "\n",
    "    def fit(self, X, y=None):\n",
    "        \"\"\"\"\"\"\n",
    "        return self\n",
    "\n",
    "    def transform(self, X, y=None):\n",
    "        return X\n",
    "\n",
    "    def fit_transform(self, X, y=None):\n",
    "        self.fit(X, y)\n",
    "        return self.transform(X, y)\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Estimator Interface\n",
    "\n",
    "Recall the estimator and predictor interfaces.\n",
    "\n",
    "```python\n",
    "class Estimator:\n",
    "\n",
    "    def fit(self, X, y=None):\n",
    "        \"\"\"Fit model to data X (and y)\"\"\"\n",
    "        self.some_attribute = self.some_fitting_method(X, y)\n",
    "        return self\n",
    "\n",
    "    def predict(self, X_test):\n",
    "        \"\"\"Make prediction based on passed features\"\"\"\n",
    "        pred = self.make_prediction(X_test)\n",
    "        return pred\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Pipelines\n",
    "\n",
    "Putting it together with a Pipeline\n",
    "\n",
    "```python\n",
    "from sklearn.pipeline import Pipeline\n",
    "\n",
    "estimator = Pipeline(steps=[\n",
    "    ('transformer1', Transformer(*args1)),\n",
    "    ('transformer2', Transformer(*args2)),\n",
    "    ('estimator', Estimator(*args))\n",
    "])\n",
    "\n",
    "estimator.fit(X_train, y_train)\n",
    "\n",
    "y_fitted = estimator.predict(X_test)\n",
    "```\n",
    "\n",
    "By chaining together transformer estimators, our code is much easier to deal with than it would have been otherwise.\n",
    "\n",
    "Under the hood, this calls fit on the first transformer, then transform on X and passes the transformed X to the next transformer until the final estimator. At the final estimator, the pipeline simply calls fit on the previously transformed X and y.\n",
    "\n",
    "Each of these estimators become a `named_step` on the Pipeline object. Available to be inspected\n",
    "\n",
    "```python\n",
    "estimator.named_steps['transformer2']\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise\n",
    "\n",
    "Take a bit and look back at what we've done in this notebook. Make a scikit-learn `Pipeline` object that uses only the comments data. Transform the comment data to tf-idf, perform TruncatedSVD, and fit a RandomForest model. What happens to the ROC-AUC when you change the number of components in the TruncatedSVD to 2? To 20? Your final classifier step should look like\n",
    "\n",
    "```python\n",
    "estimator.fit(violations_train, y_train)\n",
    "```\n",
    "\n",
    "We start here, because you're often going to need to either write your own Transformers to deal with pandas DataFrames, or you might be interested in exploring the [sklearn-pandas](https://github.com/scikit-learn-contrib/sklearn-pandas) library.\n",
    "\n",
    "You may also be interested in exploring the [FeatureUnion](http://scikit-learn.org/stable/modules/generated/sklearn.pipeline.FeatureUnion.html) transformer to write your own transformer for both the text and non-text data in the same pipeline."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Your solution here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%load solutions/sklearn_pipeline.py"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}