06_classifiers_overview-part_1.ipynb

{
  "cells": [
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# IGNORE THIS CELL WHICH CUSTOMIZES LAYOUT AND STYLING OF THE NOTEBOOK !\n",
        "%matplotlib inline\n",
        "%config InlineBackend.figure_format = 'retina'\n",
        "import warnings\n",
        "\n",
        "import matplotlib.pyplot as plt\n",
        "\n",
        "warnings.filterwarnings(\"ignore\", category=FutureWarning)\n",
        "warnings.filterwarnings(\n",
        "    \"ignore\",\n",
        "    message=\"X does not have valid feature names, but [a-zA-Z]+ was fitted with feature names\",\n",
        "    category=UserWarning,\n",
        ")\n",
        "\n",
        "warnings.filterwarnings = lambda *a, **kw: None\n",
        "from IPython.core.display import HTML\n",
        "\n",
        "HTML(open(\"custom.html\", \"r\").read())"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "# Chapter 6: An overview of classifiers, Part 1\n",
        "\n",
        "<span style=\"font-size: 150%;\">Nearest Neighbors and linear-based methods</span>"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "This script gives a quick hands-on overview of **how different types of classifiers work, their advantages and their disadvantages**. This should give you an idea of a concept behind each classifier type as well as when and which classifier type to use.\n",
        "\n",
        "For the sake of visualisation we continue with 2 dimensional data examples. For different classifiers we'll be looking at their decision surfaces. Let's start with some helper functions for that:"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "import matplotlib\n",
        "import matplotlib.pyplot as plt\n",
        "import numpy as np\n",
        "\n",
        "\n",
        "def samples_color(ilabels, colors=[\"steelblue\", \"chocolate\"]):\n",
        "    \"\"\"Return colors list from labels list given as indices.\"\"\"\n",
        "    return [colors[int(i)] for i in ilabels]\n",
        "\n",
        "\n",
        "def plot_decision_surface(\n",
        "    features_2d,\n",
        "    labels,\n",
        "    classifier,\n",
        "    preprocessing=None,\n",
        "    plt=plt,\n",
        "    marker=\".\",\n",
        "    N=100,\n",
        "    alpha=0.2,\n",
        "    colors=[\"steelblue\", \"chocolate\"],\n",
        "    title=None,\n",
        "    test_features_2d=None,\n",
        "    test_labels=None,\n",
        "    test_s=60,\n",
        "):\n",
        "    \"\"\"Plot a 2D decision surface for a already trained classifier.\"\"\"\n",
        "\n",
        "    # sanity check\n",
        "    assert len(features_2d.columns) == 2\n",
        "\n",
        "    # pandas to numpy array; get min/max values\n",
        "    xy = np.array(features_2d)\n",
        "    min_x, min_y = xy.min(axis=0)\n",
        "    max_x, max_y = xy.max(axis=0)\n",
        "\n",
        "    # create mesh of NxN points; tech: `N*1j` is spec for including max value\n",
        "    XX, YY = np.mgrid[min_x : max_x : N * 1j, min_y : max_y : N * 1j]\n",
        "    points = np.c_[XX.ravel(), YY.ravel()]  # shape: (N*N)x2\n",
        "    # points = pd.DataFrame(points, columns=[\"x\", \"y\"])\n",
        "\n",
        "    # apply scikit-learn API preprocessing\n",
        "    if preprocessing is not None:\n",
        "        points = preprocessing.transform(points)\n",
        "\n",
        "    # classify grid points\n",
        "    classes = classifier.predict(points)\n",
        "\n",
        "    # plot classes color mesh\n",
        "    ZZ = classes.reshape(XX.shape)  # shape: NxN\n",
        "    plt.pcolormesh(\n",
        "        XX,\n",
        "        YY,\n",
        "        ZZ,\n",
        "        alpha=alpha,\n",
        "        cmap=matplotlib.colors.ListedColormap(colors),\n",
        "        shading=\"auto\",\n",
        "    )\n",
        "    # plot points\n",
        "    plt.scatter(\n",
        "        xy[:, 0],\n",
        "        xy[:, 1],\n",
        "        marker=marker,\n",
        "        color=samples_color(labels, colors=colors),\n",
        "    )\n",
        "    # set title\n",
        "    if title:\n",
        "        if hasattr(plt, \"set_title\"):\n",
        "            plt.set_title(title)\n",
        "        else:\n",
        "            plt.title(title)\n",
        "    # plot test points\n",
        "    if test_features_2d is not None:\n",
        "        assert test_labels is not None\n",
        "        assert len(test_features_2d.columns) == 2\n",
        "        test_xy = np.array(test_features_2d)\n",
        "        plt.scatter(\n",
        "            test_xy[:, 0],\n",
        "            test_xy[:, 1],\n",
        "            s=test_s,\n",
        "            facecolors=\"none\",\n",
        "            linewidths=2,\n",
        "            color=samples_color(test_labels),\n",
        "        );"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Since the latest version, `sklearn` offers its own method for visualizing decision boundaries: `DecisionBoundaryDisplay`. Documentation for this method can be found here: https://scikit-learn.org/stable/modules/generated/sklearn.inspection.DecisionBoundaryDisplay.html."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Nearest Neighbors\n",
        "\n",
        "The idea is very simple: to classify a sample $x$ look for **$N$ closest samples in the training data** (by default, using the Euclidean distance) and take **majority of their labels** as a result.\n",
        "\n",
        "This method does well where the fast linear classifiers would fail, such as with the XOR dataset:"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "import pandas as pd\n",
        "\n",
        "df = pd.read_csv(\"data/xor.csv\")\n",
        "df.head(2)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "features_2d = df.loc[:, (\"x\", \"y\")]\n",
        "labelv = df[\"label\"]\n",
        "\n",
        "plt.figure(figsize=(5, 5))\n",
        "plt.scatter(features_2d.iloc[:, 0], features_2d.iloc[:, 1], color=samples_color(labelv));"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "You will find the Nearest Neighbors method in the `sklearn.neighbors` module."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "from sklearn.model_selection import train_test_split\n",
        "from sklearn.neighbors import KNeighborsClassifier\n",
        "\n",
        "X_train, X_test, y_train, y_test = train_test_split(\n",
        "    features_2d, labelv, random_state=10\n",
        ")\n",
        "X_train"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Let's use 5 neighbors to learn\n",
        "classifier = KNeighborsClassifier(n_neighbors=5)\n",
        "classifier.fit(X_train, y_train)\n",
        "\n",
        "print(\"train score: {:.2f}%\".format(100 * classifier.score(X_train, y_train)))\n",
        "print(\"test score: {:.2f}%\".format(100 * classifier.score(X_test, y_test)))\n",
        "\n",
        "plt.figure(figsize=(5, 5))\n",
        "plot_decision_surface(\n",
        "    features_2d,\n",
        "    labelv,\n",
        "    classifier,\n",
        "    test_features_2d=X_test,\n",
        "    test_labels=y_test,\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "About the plot: **the points surrounded with a circle are from the test data set** (not used for learning), all other points belong to the training data.\n",
        "\n",
        "We can query directly for the closest neighbors of a point. Let's check neighborhood of the origin:"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# distances and row indices of neighbours around (0, 0):\n",
        "dist_nn, ind_nn = classifier.kneighbors([[0, 0]])\n",
        "\n",
        "# tech: simplify dimensions\n",
        "ind_nn = ind_nn.squeeze()\n",
        "dist_nn = dist_nn.squeeze()\n",
        "\n",
        "# build data frame with neighbours\n",
        "df = X_train.iloc[ind_nn, :].copy()\n",
        "df[\"label\"] = y_train.iloc[ind_nn]\n",
        "df[\"dist\"] = dist_nn\n",
        "df"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### Exercise section\n",
        "\n",
        "Load the beers dataset and experiment with a number of neighbors (`n_neighbors`) as well as with the Manhattan distance norm `p = 1` (`2` is Euclidian distance)."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "import pandas as pd\n",
        "from sklearn.model_selection import train_test_split\n",
        "from sklearn.neighbors import KNeighborsClassifier\n",
        "from sklearn.pipeline import make_pipeline\n",
        "from sklearn.preprocessing import StandardScaler\n",
        "\n",
        "df = pd.read_csv(\"data/beers.csv\")\n",
        "features = df.iloc[:, :-1]\n",
        "labelv = df.iloc[:, -1]\n",
        "\n",
        "n_neighbors_values = [2, 3, 5, 10, 20]\n",
        "p_values = [1, 2]\n",
        "# ..."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "tags": [
          "solution"
        ]
      },
      "outputs": [],
      "source": [
        "# SOLUTION\n",
        "import pandas as pd\n",
        "from sklearn.model_selection import train_test_split\n",
        "from sklearn.neighbors import KNeighborsClassifier\n",
        "from sklearn.pipeline import make_pipeline\n",
        "from sklearn.preprocessing import StandardScaler\n",
        "\n",
        "df = pd.read_csv(\"data/beers.csv\")\n",
        "features = df.iloc[:, :-1]\n",
        "labelv = df.iloc[:, -1]\n",
        "\n",
        "\n",
        "X_train, X_test, y_train, y_test = train_test_split(features, labelv, random_state=10)\n",
        "\n",
        "p_values = [1, 2]\n",
        "n_neighbors_values = [2, 3, 5, 10, 20]\n",
        "\n",
        "for p in p_values:\n",
        "    print(f\"#### Norm L{p}\")\n",
        "    for n_neighbors in n_neighbors_values:\n",
        "\n",
        "        print(\"n_neighbors =\", n_neighbors)\n",
        "\n",
        "        pipeline = make_pipeline(\n",
        "            StandardScaler(), KNeighborsClassifier(p=p, n_neighbors=n_neighbors)\n",
        "        )\n",
        "        pipeline.fit(X_train, y_train)\n",
        "\n",
        "        print(f\"  train score: {100 * pipeline.score(X_train, y_train):.2f}%\")\n",
        "        print(f\"   test score: {100 * pipeline.score(X_test, y_test):.2f}%\")\n",
        "\n",
        "    print()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "#### Optional exercise\n",
        "\n",
        "\n",
        "Instead of choosing a number of neighbors you can also specify a radius within which samples make decision, or center of a closest class. Compare decision surface for these methods, as represented by [`RadiusNeighborsClassifier`](https://scikit-learn.org/stable/modules/generated/sklearn.neighbors.RadiusNeighborsClassifier.html#sklearn.neighbors.RadiusNeighborsClassifier) and [`NearestCentroid`](https://scikit-learn.org/stable/modules/generated/sklearn.neighbors.NearestCentroid.html#sklearn.neighbors.NearestCentroid) classifiers in the scikit-learn library.\n",
        "\n",
        "<div class=\"alert alert-block alert-info\">\n",
        "<i class=\"fa fa-info-circle\"></i>\n",
        "    Choice of an specific querying algorithm (<code>algorithm</code> parameter) becomes important with larger datasets; see: <a href=\"https://scikit-learn.org/stable/modules/neighbors.html#choice-of-nearest-neighbors-algorithm\"><em>Choice of Nearest Neighbors Algorithm</em>.</a>.\n",
        "</div>"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Few words about optimization and regularization in ML\n",
        "\n",
        "Many machine learning algorithms adapt internal weights (or similar parameters) to match predictions on the training data as good as possible.\n",
        "\n",
        "Finding these weights can be formulated as an optimzation problem which **minimizes a cost function**. Solution is usually computed in iterative improvements.\n",
        "\n",
        "<table>\n",
        "    <tr><td><img src=\"./images/cost_minimization_iterative.png\" width=400px></td></tr>\n",
        "    <tr><td><center><sub>Source: <a href=\"https://towardsdatascience.com/machine-learning-fundamentals-via-linear-regression-41a5d11f5220\">https://towardsdatascience.com/machine-learning-fundamentals-via-linear-regression-41a5d11f5220</a></sub></center></td></tr>\n",
        "</table>\n",
        "\n",
        "\n",
        "As the optimization only considers training data, taking only classification error/loss as a cost function introduces overfitting. We've seen this before in the overfitting chapter. To compensate for overfitting we can add penalty for more complex models.\n",
        "\n",
        "### The cost function\n",
        "\n",
        "<div class=\"alert alert-block alert-warning\">\n",
        "<i class=\"fa fa-info-circle\"></i> The cost function to minimize consists of an term measuring a <strong>classification loss</strong> and an additional <strong>regularization penalty</strong>:\n",
        "\n",
        "$$\\text{cost} =  \\text{classification_loss} + \\lambda \\cdot \\text{regularization_penalty}$$\n",
        "\n",
        "</div>\n",
        "\n",
        "\n",
        "The **regularization weight $\\lambda$** allows to balance out both terms and must be chosen depending on the actual algorithm and the data. In general:\n",
        "\n",
        "\n",
        "* $\\lambda$ close to `0`, $$\\text{cost} \\approx \\text{classification_loss},$$ implies more focus on training data, thus, more complex models and possible overfitting,\n",
        "\n",
        "\n",
        "* $\\lambda$ very large, $$\\text{cost} \\approx \\lambda\\cdot\\text{regularization_penalty},$$  implies less focus on training data, thus, simpler models and possible underfitting.\n",
        "\n",
        "\n",
        "Weighting-in regularization penalty relates to [Occam's razor](https://en.wikipedia.org/wiki/Occam%27s_razor) which states **_\"simpler solutions are more likely to be correct than complex ones.\"_**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Logistic Regression\n",
        "\n",
        "The name is misleading: a) despite \"logistic\" the method is linear, b) despite \"regression\" (as in unsupervised learning), it's a classification method.\n",
        "\n",
        "The method learns weights $w_1,\\cdots,w_n$ for sum of features and the threshold $b$, i.e. to learn a spearation hyper-plane:\n",
        "\n",
        "$$\n",
        "\\text{class}~0:\\quad w_1 \\cdot \\text{feature}_1 + \\ldots + w_n \\cdot \\text{feature}_n \\geq b\n",
        "$$\n",
        "$$\n",
        "\\text{class}~1:\\quad w_1 \\cdot \\text{feature}_1 + \\ldots + w_n \\cdot \\text{feature}_n \\lt b\n",
        "$$\n",
        "\n",
        "Then, to classify, transform the weighting result to a probability of belonging to one of the two classes using the **logistic function**:\n",
        "\n",
        "\n",
        "$$\n",
        "p\\left(x_1,\\cdots,x_n\\right)=\\frac{1}{1+\\exp\\left(b - \\sum_{i=1}^{n} w_i \\cdot x_i\\right)}.\n",
        "$$\n",
        "\n",
        "which looks like that:"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# classification threshold\n",
        "threshold = 3\n",
        "# results of weighted sum (around the threshold)\n",
        "weighted_sum = np.linspace(-10 + threshold, 10 + threshold, 100)\n",
        "# logistic function\n",
        "plt.plot(weighted_sum, 1 / (1 + np.exp(threshold - weighted_sum)))\n",
        "plt.axvline(x=threshold, linestyle=\"--\")\n",
        "\n",
        "plt.ylabel(\"probability that sample's class is 0\")\n",
        "plt.xlabel(\"weighted sum of features values\");"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### Demonstration\n",
        "\n",
        "Let's use a (almost) line-separable dataset:"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "import pandas as pd\n",
        "\n",
        "df = pd.read_csv(\"data/line_separable_2d.csv\")\n",
        "df.head(2)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "features_2d = df.loc[:, (\"x\", \"y\")]\n",
        "labelv = df[\"label\"]\n",
        "\n",
        "plt.figure(figsize=(5, 5))\n",
        "plt.scatter(features_2d.iloc[:, 0], features_2d.iloc[:, 1], color=samples_color(labelv));"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "You can find the Logistic Regression method in the `sklearn.linear_model` module.\n",
        "\n",
        "In scikit-learn `LogisticRegression` the regularization weight is passed here in \"inverse\", as a classification weight parameter `C` (default `1`), meaning that it multiplies the classification loss, not the regularization penalty:\n",
        "\n",
        "$$\\text{cost} =  \\text{C}\\cdot\\text{classification_loss} + \\text{regularization_penalty}$$\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "from sklearn.linear_model import LogisticRegression\n",
        "from sklearn.model_selection import train_test_split\n",
        "\n",
        "X_train, X_test, y_train, y_test = train_test_split(features_2d, labelv, random_state=0)\n",
        "\n",
        "classifier = LogisticRegression(C=1, random_state=0)\n",
        "classifier.fit(X_train, y_train)\n",
        "print(\"train score: {:.2f}%\".format(100 * classifier.score(X_train, y_train)))\n",
        "print(\"test score: {:.2f}%\".format(100 * classifier.score(X_test, y_test)))"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Feature weights are available via `.coef_` attribute, whereas threshold is the negated `.intercept_` attribute. With these we can plot separation line.\n",
        "\n",
        "Let's see how does it look like and what happens if we put more weight on the classification loss (increase `C` parameter)."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "fig, ax_arr = plt.subplots(ncols=2, nrows=1, figsize=(2 * 5, 5))\n",
        "\n",
        "plot_decision_surface(\n",
        "    features_2d,\n",
        "    labelv,\n",
        "    classifier,\n",
        "    test_features_2d=X_test,\n",
        "    test_labels=y_test,\n",
        "    plt=ax_arr[0],\n",
        "    title=\"C=1\",\n",
        ")\n",
        "\n",
        "print(\"feature weights:\", classifier.coef_)\n",
        "\n",
        "\n",
        "def plot_separation_line(features_2d, linear_classifier, plt=plt):\n",
        "    \"\"\"Plot a separation line for 2D dataset\"\"\"\n",
        "\n",
        "    assert hasattr(linear_classifier, \"coef_\")\n",
        "\n",
        "    w = linear_classifier.coef_[0]\n",
        "    b = -linear_classifier.intercept_  # NOTE: intercept = negative threshold\n",
        "\n",
        "    # separation line: w[0] * x + w[1] * y - b == 0\n",
        "    feat_x = features_2d.iloc[:, 0]\n",
        "    x = np.linspace(np.min(feat_x), np.max(feat_x), 2)\n",
        "    y = (b - w[0] * x) / w[1]\n",
        "    plt.plot(x, y, color=\"k\", linestyle=\":\")\n",
        "\n",
        "\n",
        "plot_separation_line(features_2d, classifier, plt=ax_arr[0])\n",
        "\n",
        "\n",
        "print()\n",
        "print()\n",
        "print(\"With C=100\")\n",
        "print()\n",
        "\n",
        "classifier = LogisticRegression(C=100, random_state=0)\n",
        "classifier.fit(X_train, y_train)\n",
        "print(\"train score: {:.2f}%\".format(100 * classifier.score(X_train, y_train)))\n",
        "print(\"test score: {:.2f}%\".format(100 * classifier.score(X_test, y_test)))\n",
        "print(\"feature weights:\", classifier.coef_)\n",
        "\n",
        "plot_decision_surface(\n",
        "    features_2d,\n",
        "    labelv,\n",
        "    classifier,\n",
        "    test_features_2d=X_test,\n",
        "    test_labels=y_test,\n",
        "    plt=ax_arr[1],\n",
        "    title=\"C=100\",\n",
        ")\n",
        "plot_separation_line(features_2d, classifier, plt=ax_arr[1])"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### Exercise section\n",
        "\n",
        "1. Why does the test score drop when we penalize more misclassifications?\n",
        "2. For the higher dimensional beers dataset experiment with both `C` and `penalty` parameters of the linear regression classfier. Compare scores and the resulting weights. What does the `l1` penalty do? What is the sweet spot of the \"inverse regularization\" `C`?\n",
        "  "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "tags": [
          "solution"
        ]
      },
      "outputs": [],
      "source": [
        "# SOLUTION 1\n",
        "\n",
        "# With C=100 we try hard to get all training points correctly classified, whereas with C=1\n",
        "# you can see we allow misclassification in training, in order to possibly get more general\n",
        "# model and avoid overfitting.\n",
        "#\n",
        "# You can see in the test data, that reverse - one misclassfied point with C=100.\n",
        "# If we would have that point for for training, the line would look more like in C=1 case.\n",
        "# (Go regularization!)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "import pandas as pd\n",
        "from sklearn.linear_model import LogisticRegression\n",
        "from sklearn.model_selection import train_test_split\n",
        "from sklearn.pipeline import make_pipeline\n",
        "from sklearn.preprocessing import StandardScaler\n",
        "\n",
        "df = pd.read_csv(\"data/beers.csv\")\n",
        "features = df.iloc[:, :-1]\n",
        "labelv = df.iloc[:, -1]\n",
        "\n",
        "C_values = [0.01, 0.1, 1, 10, 100, 1000]\n",
        "penalty_values = [\"l1\", \"l2\"]\n",
        "# ..."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "tags": [
          "solution"
        ]
      },
      "outputs": [],
      "source": [
        "# SOLUTION 2\n",
        "import pandas as pd\n",
        "from sklearn.linear_model import LogisticRegression\n",
        "from sklearn.model_selection import train_test_split\n",
        "from sklearn.pipeline import make_pipeline\n",
        "from sklearn.preprocessing import StandardScaler\n",
        "\n",
        "df = pd.read_csv(\"data/beers.csv\")\n",
        "print(df.head(2))\n",
        "\n",
        "features = df.iloc[:, :-1]\n",
        "labelv = df.iloc[:, -1]\n",
        "\n",
        "X_train, X_test, y_train, y_test = train_test_split(features, labelv, random_state=10)\n",
        "\n",
        "C_values = [0.01, 0.1, 1, 10, 100, 1000]\n",
        "penalty_values = [\"l1\", \"l2\"]\n",
        "\n",
        "print()\n",
        "for norm in penalty_values:\n",
        "    print(\"#### Norm\", norm)\n",
        "    print()\n",
        "    for C in C_values:\n",
        "        print(\"C:\", C)\n",
        "        # Note: use non-default solver for L1 penalty support\n",
        "        # Note: increase max iterations 10x for solver's convergence\n",
        "        pipeline = make_pipeline(\n",
        "            StandardScaler(),\n",
        "            LogisticRegression(\n",
        "                C=C, solver=\"liblinear\", penalty=norm, dual=False, max_iter=10000\n",
        "            ),\n",
        "        )\n",
        "        pipeline.fit(X_train, y_train)\n",
        "        print(f\"  train score: {100 * pipeline.score(X_train, y_train):.2f}%\")\n",
        "        print(f\"   test score: {100 * pipeline.score(X_test, y_test):.2f}%\")\n",
        "        print(\"      weights:\", pipeline[-1].coef_[0])\n",
        "    print()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "<div class=\"alert alert-block alert-info\">\n",
        "\n",
        "<p><i class=\"fa fa-info-circle\"></i>&nbsp;\n",
        "The <strong>classification loss</strong> in logistic regression is a so called <em>negative-log likelihood</em>, i.e. a negative logarithm of the logistic probability above:\n",
        "<p/>\n",
        "    \n",
        "<p>\n",
        "$$ \\text{classification_loss} = -\\log(p(x^k; p^k)) = \\log{\\left(1+\\exp{\\left(y^k\\left(b - \\sum_{i=1}^{n} w_i x_i^k\\right)\\right)}\\right)}$$\n",
        "<p/>\n",
        "\n",
        "<p>\n",
        "where $y^k$ is -1 or 1, representing class of $k$-th sample from the training data, corresponding, respectively, to class below and above the threshold (the separation line).\n",
        "\n",
        "The $+/-$ sign for the class penalizes missclassifications. If sample is below the threshold $\\sum_{i=1}^{n} w_i x_i^k < b$ and have the correct class $y^k = -1$, then we have $\\exp{\\left(\\text{negative value}\\right)}$ giving small loss. In case of misclassification $\\exp{\\left(\\text{positive value}\\right)}$ gives a much bigger loss.\n",
        "</p>\n",
        "</div>\n",
        "\n",
        "<div class=\"alert alert-block alert-info\">\n",
        "<p><i class=\"fa fa-info-circle\"></i>&nbsp;\n",
        "The <strong>reqularization penalty</strong> in logistic regression is a <em>norm of the learnt weights</em>, denoted as:\n",
        "\n",
        "<p>\n",
        "$$\\text{regularization_penalty} = \\left\\lVert w \\right\\rVert_p$$\n",
        "</p>\n",
        "\n",
        "<p>\n",
        "Using <em>L1 norm</em> ($p=1$, Manhattan distance from origin, which is sum of absolute weight values) is know for finding sparse solutions, i.e. eliminating features (weight equal to 0) when they are have low significance. With the default <em>L2 norm</em> ($p=2$, Euclidian distance from origin, which is square root of sum of squared weight values), weights of insignificant features would have small non-zero values instead.\n",
        "</p>\n",
        "\n",
        "<p>\n",
        "In <code>LogisticRegression</code> class, <code>penalty</code> parameter allows to specify type of norm to use.\n",
        "</p>\n",
        "\n",
        "<p>\n",
        "Note that any solution weights and its threshold can be scaled to give the same result. Thus the regularization penalty not only prevents overfitting but also ensures a unique solution.\n",
        "</p>\n",
        "\n",
        "</div>"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Linear SVM\n",
        "\n",
        "Support-Vector Machine (SVM) classifier tries to separate two classes with a line by **finding data points (support vectors) lying closest to the separation plane**. These points determine separation plane (weights and threshold/intercept).\n",
        "\n",
        "The weights are learned such that the **margin between support vectors of different classes is maximized**.\n",
        "\n",
        "<table>\n",
        "    <tr><td><img src=\"./images/svm_margin.png\" width=400px></td></tr>\n",
        "    <tr><td><center><sub>Source: <a href=\"https://en.wikipedia.org/wiki/Support-vector_machine\">https://en.wikipedia.org/wiki/Support-vector_machine</a></sub></center></td></tr>\n",
        "</table>\n",
        "\n",
        "Like in linear regression the classification is based on a weighted sum of the features (and margin maximization corresponds to minimization of the regularization penalty). \n",
        "\n",
        "Analogously to the Nearest Neighbors method the data points (support vectors) decide the class of a new data sample."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### Demonstration\n",
        "\n",
        "Let's try it out on the line-separable dataset.\n",
        "\n",
        "You will find the SVM method in the `sklearn.svm` module."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "import pandas as pd\n",
        "from sklearn.model_selection import train_test_split\n",
        "from sklearn.svm import LinearSVC\n",
        "\n",
        "df = pd.read_csv(\"data/line_separable_2d.csv\")\n",
        "features_2d = df.loc[:, (\"x\", \"y\")]\n",
        "labelv = df[\"label\"]\n",
        "\n",
        "X_train, X_test, y_train, y_test = train_test_split(features_2d, labelv, random_state=0)\n",
        "\n",
        "classifier = LinearSVC(C=1)\n",
        "classifier.fit(X_train, y_train)\n",
        "print(\"train score: {:.2f}%\".format(100 * classifier.score(X_train, y_train)))\n",
        "print(\"test score: {:.2f}%\".format(100 * classifier.score(X_test, y_test)))"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Again, let's see how does the separation line look like here and what happens if we put more weight on the classification loss (increase `C` parameter)."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "fig, ax_arr = plt.subplots(ncols=2, nrows=1, figsize=(2 * 5, 5))\n",
        "\n",
        "plot_decision_surface(\n",
        "    features_2d,\n",
        "    labelv,\n",
        "    classifier,\n",
        "    test_features_2d=X_test,\n",
        "    test_labels=y_test,\n",
        "    plt=ax_arr[0],\n",
        "    title=\"C=1\",\n",
        ")\n",
        "\n",
        "print(\"feature weights:\", classifier.coef_)\n",
        "\n",
        "\n",
        "def plot_margins(features_2d, linear_classifier, plt=plt):\n",
        "    \"\"\"Plot a separation line and margin lines for 2D dataset\"\"\"\n",
        "\n",
        "    assert hasattr(linear_classifier, \"coef_\")\n",
        "\n",
        "    w = linear_classifier.coef_[0]\n",
        "    b = -linear_classifier.intercept_  # NOTE: intercept = negative threshold\n",
        "\n",
        "    # separation line: w[0] * x + w[1] * y - b == 0\n",
        "    feat_x = features_2d.iloc[:, 0]\n",
        "    x = np.linspace(np.min(feat_x), np.max(feat_x), 2)\n",
        "    y = (b - w[0] * x) / w[1]\n",
        "    plt.plot(x, y, color=\"k\", linestyle=\":\")\n",
        "\n",
        "    # margin lines: w[0] * x + w[1] * y - b == +/-1\n",
        "    y = ((b - 1) - w[0] * x) / w[1]\n",
        "    plt.plot(x, y, color=\"r\", linestyle=\":\")\n",
        "    y = ((b + 1) - w[0] * x) / w[1]\n",
        "    plt.plot(x, y, color=\"r\", linestyle=\":\")\n",
        "\n",
        "\n",
        "plot_margins(features_2d, classifier, plt=ax_arr[0])\n",
        "\n",
        "\n",
        "print()\n",
        "print()\n",
        "print(\"With C=100\")\n",
        "print()\n",
        "\n",
        "# higher C = more narrow (\"harder\") margin\n",
        "# Note: increase max iterations 50x for solver's convergence\n",
        "classifier = LinearSVC(C=100, max_iter=50000)\n",
        "classifier.fit(X_train, y_train)\n",
        "print(\"train score: {:.2f}%\".format(100 * classifier.score(X_train, y_train)))\n",
        "print(\"test score: {:.2f}%\".format(100 * classifier.score(X_test, y_test)))\n",
        "print(\"feature weights:\", classifier.coef_)\n",
        "\n",
        "plot_decision_surface(\n",
        "    features_2d,\n",
        "    labelv,\n",
        "    classifier,\n",
        "    test_features_2d=X_test,\n",
        "    test_labels=y_test,\n",
        "    plt=ax_arr[1],\n",
        "    title=\"C=100\",\n",
        ")\n",
        "plot_margins(features_2d, classifier, plt=ax_arr[1]);"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Why are there **training data samples within the margin**?\n",
        "\n",
        "This is because the scikit-learn `LinearSVC` implementation by default uses so called **\"soft margins\"**, (\"hinge\" loss function in `loss` parameter). They allow support vectors to lie within the -1, 1 margin (with appropriately lower weights than -1, 1).\n",
        "\n",
        "You can control \"softness\" or \"hardness\" of classification loss by, respectively, decreasing or increasing its weight (parameter `C` of the `LinearSVC` class)."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### Exercise section\n",
        "\n",
        "1. It looks like we did train our classifier \"perfectly\" with \"harder\" margins. Why is the score then lower then previously?\n",
        "2. For the higher dimensional beers dataset experiment with both `C` and `penalty` parameters of the linear SVM classfier (note: set `dual=False` to work with `penalty='l1'`). Compare scores and the resulting weights.\n",
        "  "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "import pandas as pd\n",
        "from sklearn.model_selection import train_test_split\n",
        "from sklearn.pipeline import make_pipeline\n",
        "from sklearn.preprocessing import StandardScaler\n",
        "from sklearn.svm import LinearSVC\n",
        "\n",
        "df = pd.read_csv(\"data/beers.csv\")\n",
        "features = df.iloc[:, :-1]\n",
        "labelv = df.iloc[:, -1]\n",
        "\n",
        "C_values = [0.01, 0.1, 1, 10, 100, 1000]\n",
        "penalty_values = [\"l1\", \"l2\"]\n",
        "\n",
        "# ..."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "tags": [
          "solution"
        ]
      },
      "outputs": [],
      "source": [
        "# SOLUTION 1\n",
        "\n",
        "# Again, with C=100 we've just tried to hard to get all training points correctly classified,\n",
        "# but this time it meant essentially no points within the margin. Thus, by overfitting we\n",
        "# lost the linear trend in the data, which is represented by the one test data sample that\n",
        "# just did not make it over the separation line (and a bit overrepresented by the other quite\n",
        "# badly misclassfied test sample)."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "tags": [
          "solution"
        ]
      },
      "outputs": [],
      "source": [
        "# SOLUTION 2\n",
        "import pandas as pd\n",
        "from sklearn.model_selection import train_test_split\n",
        "from sklearn.pipeline import make_pipeline\n",
        "from sklearn.preprocessing import StandardScaler\n",
        "from sklearn.svm import LinearSVC\n",
        "\n",
        "df = pd.read_csv(\"data/beers.csv\")\n",
        "print(df.head(2))\n",
        "\n",
        "features = df.iloc[:, :-1]\n",
        "labelv = df.iloc[:, -1]\n",
        "\n",
        "X_train, X_test, y_train, y_test = train_test_split(features, labelv, random_state=10)\n",
        "\n",
        "C_values = [0.01, 0.1, 1, 10, 100, 1000]\n",
        "penalty_values = [\"l1\", \"l2\"]\n",
        "\n",
        "print()\n",
        "for norm in penalty_values:\n",
        "    print(\"#### Norm\", norm)\n",
        "    print()\n",
        "    for C in C_values:\n",
        "        print(\"C:\", C)\n",
        "        # Note: increase max iterations 10x for solver's convergence\n",
        "        pipeline = make_pipeline(\n",
        "            StandardScaler(),\n",
        "            LinearSVC(C=C, penalty=norm, dual=False, max_iter=10000),\n",
        "        )\n",
        "        pipeline.fit(X_train, y_train)\n",
        "        print(f\"  train score: {100 * pipeline.score(X_train, y_train):.2f}%\")\n",
        "        print(f\"   test score: {100 * pipeline.score(X_test, y_test):.2f}%\")\n",
        "        print(\"      weights:\", pipeline[-1].coef_[0])\n",
        "    print()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Kernel based SVM\n",
        "\n",
        "### Dealing with not linearily-separable data\n",
        "\n",
        "Data is usually not at all linearily separable.\n",
        "\n",
        "Let's revisit the \"circle\" data."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "import pandas as pd\n",
        "\n",
        "df = pd.read_csv(\"data/circle.csv\")\n",
        "df.head(2)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "features_2d = df.loc[:, (\"x\", \"y\")]\n",
        "labelv = df[\"label\"]\n",
        "\n",
        "plt.figure(figsize=(5, 5))\n",
        "plt.scatter(features_2d.iloc[:, 0], features_2d.iloc[:, 1], color=samples_color(labelv));"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Linear methods like SVM fail in such cases to learn the underlying pattern:"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "from sklearn.model_selection import train_test_split\n",
        "from sklearn.svm import LinearSVC\n",
        "\n",