14_LogisticRegression.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "executionInfo": {
     "elapsed": 863,
     "status": "ok",
     "timestamp": 1605029011180,
     "user": {
      "displayName": "Adam Gilbert",
      "photoUrl": "https://lh3.googleusercontent.com/a-/AOh14Gj1KWqLKG2SAQxx8KvqzqhM94EvqxkXs2iOw7gQWQ=s64",
      "userId": "02074648007699947871"
     },
     "user_tz": 300
    },
    "id": "7bXwjMhVMeea"
   },
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "\n",
    "from sklearn import datasets\n",
    "from sklearn.model_selection import train_test_split\n",
    "from sklearn.linear_model import LogisticRegression\n",
    "from sklearn.model_selection import cross_val_score\n",
    "from sklearn.metrics import confusion_matrix\n",
    "\n",
    "churn = pd.read_csv(\"https://raw.githubusercontent.com/rajeevratan84/datascienceforbusiness/master/WA_Fn-UseC_-Telco-Customer-Churn.csv\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "YtrVKOOVNzii"
   },
   "source": [
    " # Logistic regression for binary classification\n",
    "\n",
    "Last week we explored the utility of applying decision trees as classification models. We applied decision tree classifiers in two cases, (i) to predict a binary response variable indicating whether or not a customer at a subscription-based firm would *churn* or not, and then (ii) a computer vision problem to classify handwritten digits -- there were 10 different classes in this second example. We saw that decision tree classifiers worked fine in both the binary and multiclass cases without changing any of our workflow. We did note, however, that a drawback to decision tree classifiers that we did not necessarily obtain a *propensity* estimate (an estimate of the probability with which we believe an observation belongs to the class it was assigned to). In this week's notebook we explore *logistic regression*, a classification model which does indeed provide us with a propensity estimate but which has the drawback of only being applicable in binary classification. We'll revisit the customer *churn* example here.\n",
    "\n",
    "At the beginning of last week's notebook on decision tree classifiers we noted that a linear regression model could not be applied to a classification problem. All linear regression models aside from the null model ($\\mathbb{E}\\left[y\\right] = \\beta_0$) become unbounded. Logistic regression is indeed a regression model -- it predicts a numerical output, but that numerical output is bounded between 0 and 1.\n",
    "\n",
    "This notebook introduces logistic regression. While logistic regression is not well-suited to the multi-class problem, it is a really common method applied to the two-class case. Logistic regression even provides some potential advantages over a more flexible model, like a tree. The output of a logistic model is a propensity for belonging to the class corresponding to an output of 1 -- the form for a logistic regression model is expressed below.\n",
    "$$\\mathbb{P}\\left[Y = 1\\right] = \\displaystyle{\\frac{e^{\\beta_0 +\\beta_1 x_1 + \\beta_2 x_2 + ... + \\beta_k x_k}}{1 + e^{\\beta_0 +\\beta_1 x_1 + \\beta_2 x_2 + ... + \\beta_k x_k}}}$$\n",
    "\n",
    "Notice that the exponentials are being \"raised to [one of our usual] regression models\".\n",
    "\n",
    "## A reminder on the `churn` data\n",
    "\n",
    "We are bringing back last week's first classification dataset which dealt with customer *churn*. As a reminder, churn is when an existing customer, user, player, subscriber or any kind of return client stops doing business or ends the relationship with a company. A natural goal is to figure our which customers may be most likely churn in future so that we can provide preventative outreach.\n",
    "\n",
    "We'll start by preparing our data in the same way it was prepared for last week's decision tree classifier."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 258
    },
    "executionInfo": {
     "elapsed": 648,
     "status": "ok",
     "timestamp": 1605029580560,
     "user": {
      "displayName": "Adam Gilbert",
      "photoUrl": "https://lh3.googleusercontent.com/a-/AOh14Gj1KWqLKG2SAQxx8KvqzqhM94EvqxkXs2iOw7gQWQ=s64",
      "userId": "02074648007699947871"
     },
     "user_tz": 300
    },
    "id": "SlzvtgN5Meem",
    "outputId": "7fb7522b-f569-494c-d348-5c16088f78cc"
   },
   "outputs": [],
   "source": [
    "churn = churn[churn[\"TotalCharges\"] != \" \"]\n",
    "churn[\"TotalCharges\"] = pd.to_numeric(churn[\"TotalCharges\"])\n",
    "\n",
    "num_cols = [\"tenure\", \"MonthlyCharges\", \"TotalCharges\"]\n",
    "unique_id = [\"customerID\"]\n",
    "cat_cols = [name for name in list(churn.columns) if ((name not in num_cols) & (name not in unique_id))]\n",
    "\n",
    "churn = pd.get_dummies(churn, columns = cat_cols, drop_first = True)\n",
    "\n",
    "churn.drop([\"customerID\"], axis = 1, inplace = True)\n",
    "\n",
    "X = churn.drop([\"Churn_Yes\"], axis = 1)\n",
    "y = churn[\"Churn_Yes\"]\n",
    "\n",
    "X_train, X_safe, y_train, y_safe = train_test_split(X, y, test_size = 0.1, random_state = 42)\n",
    "\n",
    "X_train.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "kNatwzn4Meeq"
   },
   "source": [
    "## Logistic Regression\n",
    "\n",
    "As a recap, a logistic regression model can be used to predict a two-class response variable (for example: *churn* or *no churn*). The outputs of a logistic regression model are bounded between 0 and 1, and their outputs can be interpreted as the propensity (likelihood) for an observation to belong to the class assigned by 1. The form for a logistic regression model is expressed below.\n",
    "$$\\mathbb{P}\\left[Y = 1\\right] = \\displaystyle{\\frac{e^{\\beta_0 +\\beta_1 x_1 + \\beta_2 x_2 + ... + \\beta_k x_k}}{1 + e^{\\beta_0 +\\beta_1 x_1 + \\beta_2 x_2 + ... + \\beta_k x_k}}}$$\n",
    "\n",
    "Notice that the exponentials are being \"raised to [one of our usual] regression models\".\n",
    "\n",
    "## Building a Logistic Regression Model\n",
    "\n",
    "Now that our data has been pre-processed and split, we'll fit our first logistic regression model."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "executionInfo": {
     "elapsed": 925,
     "status": "ok",
     "timestamp": 1605029640160,
     "user": {
      "displayName": "Adam Gilbert",
      "photoUrl": "https://lh3.googleusercontent.com/a-/AOh14Gj1KWqLKG2SAQxx8KvqzqhM94EvqxkXs2iOw7gQWQ=s64",
      "userId": "02074648007699947871"
     },
     "user_tz": 300
    },
    "id": "yb9Gmsla5VVY",
    "outputId": "e22395dd-2a44-4ca1-d66c-cdcb02d083e5"
   },
   "outputs": [],
   "source": [
    "#we define the model object\n",
    "log_reg_clf = LogisticRegression()\n",
    "#and now we fit the model\n",
    "log_reg_clf.fit(X_train, y_train)\n",
    "\n",
    "print(\"Intercept: \", log_reg_clf.intercept_[0])\n",
    "print(pd.DataFrame({\"Variable\" : X_train.columns, \"Coefficient\" : log_reg_clf.coef_[0]}))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "GS6UuXSeMeex"
   },
   "source": [
    "From the model output above we can see that the following characteristics **increase** the likelihood of *churn* (because their coefficients are positive):\n",
    "+ The higher the customers `MonthlyCharges` and `TotalCharges`, the more likely they are to *churn*.\n",
    "+ If they are a senior citizen.\n",
    "+ The customer has a partner.\n",
    "+ If they have no phone line (or if they have multiple phone lines).\n",
    "+ If they have fiber-optic internet.\n",
    "+ If they stream TV and/or Movies.\n",
    "+ If they have paperless billing.\n",
    "+ If they pay their bill by electronic check.\n",
    "\n",
    "If we decide that our model is a good one, we can look for profiles fitting these characteristics to prioritize outreach. Check back to last week's notebook -- do those decision trees indicate the same sorts of trends?\n",
    "\n",
    "We'll continue now with evaluating our model. Similar to last week, we'll use a confusion matrix along with the *accuracy*, *precision*, and *recall* metrics to assess our model utility. We don't typically do this, but we will utilize the training data for computing these metrics since this is just our baseline model. We'll use cross-validation shortly in order to tune hyperparameters and compare a few different logistic regression models."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "executionInfo": {
     "elapsed": 593,
     "status": "ok",
     "timestamp": 1605030625480,
     "user": {
      "displayName": "Adam Gilbert",
      "photoUrl": "https://lh3.googleusercontent.com/a-/AOh14Gj1KWqLKG2SAQxx8KvqzqhM94EvqxkXs2iOw7gQWQ=s64",
      "userId": "02074648007699947871"
     },
     "user_tz": 300
    },
    "id": "5AYyJovGMeey",
    "outputId": "5cb2e785-ad18-47a2-bd38-b53cc56b5716"
   },
   "outputs": [],
   "source": [
    "confusion_matrix(y_train, log_reg_clf.predict(X_train))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "dOdntau-Mee4"
   },
   "source": [
    "From the confusion matrix, we see that our accuracy is about 80.69% while our recall is 55.01% and precision is 66.33%. Remember that our goal is not necessarily to build a highly accurate model. In this particular application, we really want a model with high recall so that we can outreach customers who are at risk of *churn* before we lose their business. These outreaches will likely come along with promotional offers, so while we want a high recall metric, we don't want to lose too much precision because that will lead to lost revenue. We have a challenge -- **increase recall without decreasing precision** [too much]!\n",
    "\n",
    "## Competing models\n",
    "\n",
    "We'll build a set of competing models here and try using cross-validation to help us identify a model which will have high recall without leading to poor precision.\n",
    "\n",
    "### Hyperparameters for logistic regression\n",
    "\n",
    "Logistic regression has several hyperparameters which can be set. I'll list some of the ones we'll use below, but you can [see an exhaustive list here](https://scikit-learn.org/stable/modules/generated/sklearn.linear_model.LogisticRegression.html). \n",
    "+ The `C` hyperparameter essentially supplies a budget which can be spent on coefficients.\n",
    "+ The `penalty` hyperparameter determines how we compute our spent budget. The details are beyond the scope of our course but you can look up `l1` and `l2` regularization, or take MAT300 to learn more.\n",
    "+ Both the `C` and `penalty` hyperparameters can be used to reduce overfitting."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "executionInfo": {
     "elapsed": 11792,
     "status": "ok",
     "timestamp": 1605031343823,
     "user": {
      "displayName": "Adam Gilbert",
      "photoUrl": "https://lh3.googleusercontent.com/a-/AOh14Gj1KWqLKG2SAQxx8KvqzqhM94EvqxkXs2iOw7gQWQ=s64",
      "userId": "02074648007699947871"
     },
     "user_tz": 300
    },
    "id": "jQCRBFzkMee5",
    "outputId": "26e016b5-c6c9-45ea-f5d9-e20d116b942e"
   },
   "outputs": [],
   "source": [
    "budget = [0.1, 0.5, 1, 10, 100, 1000]\n",
    "regularization = [\"l1\", \"l2\"]\n",
    "\n",
    "for budget in budget:\n",
    "    for reg in regularization:\n",
    "        log_reg_clf = LogisticRegression(C = budget, penalty = reg, solver = \"liblinear\")\n",
    "        \n",
    "        cv_scores = cross_val_score(log_reg_clf, X_train, y_train, cv = 10, scoring = \"recall\")\n",
    "        \n",
    "        print(\"Budget \", budget, \" and \", reg, \" regularization:\")\n",
    "        print(\"\\t recall rates: \", cv_scores)\n",
    "        print(\"\\t cv-average-recall: \", cv_scores.mean())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "mFY1b6DxMee8"
   },
   "source": [
    "From the regression output above, we see that `l1` regularization is superior to `l2` regularization here but that recall continues to improve with larger and larger budgets. Let's run another round of cross-validation, but we'll include only `l1` regularization and allow for larger coefficient budgets."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "executionInfo": {
     "elapsed": 8769,
     "status": "ok",
     "timestamp": 1605031600841,
     "user": {
      "displayName": "Adam Gilbert",
      "photoUrl": "https://lh3.googleusercontent.com/a-/AOh14Gj1KWqLKG2SAQxx8KvqzqhM94EvqxkXs2iOw7gQWQ=s64",
      "userId": "02074648007699947871"
     },
     "user_tz": 300
    },
    "id": "tfoy3LUHMee9",
    "outputId": "0d19cba5-760d-4a61-bea7-9cb58ca025df"
   },
   "outputs": [],
   "source": [
    "largerBudgets = [10**3, 10**4, 10**5, 10**6, 10**7, 10**8]\n",
    "\n",
    "for budget in largerBudgets:\n",
    "    log_reg_clf = LogisticRegression(C = budget, penalty = \"l1\", solver = \"liblinear\")\n",
    "    \n",
    "    cv_scores = cross_val_score(log_reg_clf, X_train, y_train, cv = 10, scoring = \"recall\")\n",
    "    print(\"Budget: \", budget)\n",
    "    print(\"\\t recall rates: \", cv_scores)\n",
    "    print(\"\\t cv-average-recall: \", cv_scores.mean())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "2xjctExEMefC"
   },
   "source": [
    "Okay, according to our recall rates, it looks like a budget of around `C = 10000` produces the best average recall. Since this is the case, let's re-fit that model and analyse our results."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "executionInfo": {
     "elapsed": 558,
     "status": "ok",
     "timestamp": 1605031638269,
     "user": {
      "displayName": "Adam Gilbert",
      "photoUrl": "https://lh3.googleusercontent.com/a-/AOh14Gj1KWqLKG2SAQxx8KvqzqhM94EvqxkXs2iOw7gQWQ=s64",
      "userId": "02074648007699947871"
     },
     "user_tz": 300
    },
    "id": "rVjMqyMcMefE",
    "outputId": "9ef6462d-cb6c-4cb9-c11e-92d2e31a6865"
   },
   "outputs": [],
   "source": [
    "log_reg_clf = LogisticRegression(C = 1000, penalty = \"l1\", solver = \"liblinear\")\n",
    "\n",
    "log_reg_clf.fit(X_train, y_train)\n",
    "\n",
    "confusion_matrix(y_safe, log_reg_clf.predict(X_safe))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "executionInfo": {
     "elapsed": 440,
     "status": "ok",
     "timestamp": 1605031737744,
     "user": {
      "displayName": "Adam Gilbert",
      "photoUrl": "https://lh3.googleusercontent.com/a-/AOh14Gj1KWqLKG2SAQxx8KvqzqhM94EvqxkXs2iOw7gQWQ=s64",
      "userId": "02074648007699947871"
     },
     "user_tz": 300
    },
    "id": "tJiyJ62FMefI",
    "outputId": "af00a97f-a93d-4de2-cec9-2fe5e16e1c06"
   },
   "outputs": [],
   "source": [
    "98/(54 + 98)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "bx6ruzg8MefM"
   },
   "source": [
    "Okay, on the unseen *safe* data, our final model had about 78.84% accuracy, a recall of about 50.78%, and a precision of about 64.47%. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "e7nRJfZ5Q9xU"
   },
   "source": [
    "## Closing for Logistic Regression\n",
    "\n",
    "So there it is, you've constructed and evaluated a series of logistic regression models, tuning hyperparameters along the way. There's much more to learn about logistic regression and I really encourage you to read about this classifier in Chapter 10 of our textbook.\n",
    "\n",
    "Before we leave this topic and move to uplift and ensemble methods, there are two major takeaways you should consider here.\n",
    "\n",
    "+ **Weakness:** A major weakness of logistic regression models is that they are not well-suited to classification scenarios with more than two classes. There are methods for applying them to these cases though -- the easiest method is to decide that we only really care about whether an observation belongs to a single class of interest or not -- from there we have a *one-versus-all* classification problem, which is binary. If we do indeed care about differentiating all classes from one another we could engage in building a series of *one-versus-one* classifiers and then using them in a pseudo-ensemble fashion. The challenge here, however, is that even with just four classes we end up needing to work with $4\\left(3\\right)/2 = 6$ individual models.\n",
    "+ **Strength:** Rather than outputting a *class prediction*, a logistic regression model outputs a likelihood that an observation belongs to the class represented by 1. This means that it is easy to adjust the threshold for belonging to that class -- is greater than 50% likelihood the right choice? Or would we be better off using a really strict threshold like 90%, or a really loose threshold like 20%?$^\\dagger$ Furthermore, since we are predicting a likelihood of belonging to a particular class, it is really easy to identify the predictions we are \"most certain\" about.\n",
    "\n",
    "<center><span style = \"font-size:8pt\"><i>$^\\dagger$How do you choose the appropriate threshold??? -- the threshold is a hyperparameter -- keep calm and cross-validate!</i></span></center>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "executionInfo": {
     "elapsed": 602,
     "status": "ok",
     "timestamp": 1605031858836,
     "user": {
      "displayName": "Adam Gilbert",
      "photoUrl": "https://lh3.googleusercontent.com/a-/AOh14Gj1KWqLKG2SAQxx8KvqzqhM94EvqxkXs2iOw7gQWQ=s64",
      "userId": "02074648007699947871"
     },
     "user_tz": 300
    },
    "id": "YRMzD4aDDlZa",
    "outputId": "15f0abc4-6e90-4b0e-bb9f-a084369e844c"
   },
   "outputs": [],
   "source": [
    "log_reg_clf.predict_proba(X_safe)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "9WL_O2qlMefN"
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "colab": {
   "name": "13_LogisticRegression.ipynb",
   "provenance": []
  },
  "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.8.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}