labs/linear-models/regularization-logi.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Exploring regularization for logistic regression"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Goal\n",
    "\n",
    "The goal of this lab is to explore the effect of regularization on the coefficients and accuracy of logistic regression models for a toy (wine) dataset"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Set up"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "\n",
    "np.random.seed(999)\n",
    "\n",
    "from sklearn.linear_model import LinearRegression, LogisticRegression, Lasso, Ridge\n",
    "from sklearn.ensemble import RandomForestRegressor\n",
    "from sklearn.datasets import load_wine\n",
    "from sklearn.model_selection import train_test_split\n",
    "from sklearn.metrics import confusion_matrix, precision_score, recall_score, log_loss, mean_absolute_error\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib as mpl\n",
    "%config InlineBackend.figure_format = 'retina'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "from pandas.api.types import is_numeric_dtype\n",
    "def normalize(X): \n",
    "    for colname in X.columns:\n",
    "        if is_numeric_dtype(X[colname]):\n",
    "            u = np.mean(X[colname])\n",
    "            s = np.std(X[colname])\n",
    "            X[colname] = (X[colname] - u) / s"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Load data, create 2-class problem"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "130 records for classes {0,1} from 178 records\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>alcohol</th>\n",
       "      <th>malic_acid</th>\n",
       "      <th>ash</th>\n",
       "      <th>alcalinity_of_ash</th>\n",
       "      <th>magnesium</th>\n",
       "      <th>total_phenols</th>\n",
       "      <th>flavanoids</th>\n",
       "      <th>nonflavanoid_phenols</th>\n",
       "      <th>proanthocyanins</th>\n",
       "      <th>color_intensity</th>\n",
       "      <th>hue</th>\n",
       "      <th>od280/od315_of_diluted_wines</th>\n",
       "      <th>proline</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>14.23</td>\n",
       "      <td>1.71</td>\n",
       "      <td>2.43</td>\n",
       "      <td>15.6</td>\n",
       "      <td>127.0</td>\n",
       "      <td>2.80</td>\n",
       "      <td>3.06</td>\n",
       "      <td>0.28</td>\n",
       "      <td>2.29</td>\n",
       "      <td>5.64</td>\n",
       "      <td>1.04</td>\n",
       "      <td>3.92</td>\n",
       "      <td>1065.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>13.20</td>\n",
       "      <td>1.78</td>\n",
       "      <td>2.14</td>\n",
       "      <td>11.2</td>\n",
       "      <td>100.0</td>\n",
       "      <td>2.65</td>\n",
       "      <td>2.76</td>\n",
       "      <td>0.26</td>\n",
       "      <td>1.28</td>\n",
       "      <td>4.38</td>\n",
       "      <td>1.05</td>\n",
       "      <td>3.40</td>\n",
       "      <td>1050.0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   alcohol  malic_acid   ash  alcalinity_of_ash  magnesium  total_phenols  \\\n",
       "0    14.23        1.71  2.43               15.6      127.0           2.80   \n",
       "1    13.20        1.78  2.14               11.2      100.0           2.65   \n",
       "\n",
       "   flavanoids  nonflavanoid_phenols  proanthocyanins  color_intensity   hue  \\\n",
       "0        3.06                  0.28             2.29             5.64  1.04   \n",
       "1        2.76                  0.26             1.28             4.38  1.05   \n",
       "\n",
       "   od280/od315_of_diluted_wines  proline  \n",
       "0                          3.92   1065.0  \n",
       "1                          3.40   1050.0  "
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "wine = load_wine()\n",
    "df_wine = pd.DataFrame(data=wine.data, columns=wine.feature_names)\n",
    "df_wine['y'] = wine.target\n",
    "df_wine = df_wine[df_wine['y']<2] # Only do 2-class problem {0,1}\n",
    "\n",
    "X = df_wine.drop('y', axis=1)\n",
    "y = df_wine['y']\n",
    "print(f\"{len(X)} records for classes {{0,1}} from {len(wine.data)} records\")\n",
    "X.head(2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Examine coefficients of nonregularized OLS model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**1. Split off validation/test set and train unregularized linear regression model**\n",
    "\n",
    "Select a seed and random state known to yield poor validation set accuracy (depending on the split, scores will fluctuate, particularly in the presence of outliers.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, shuffle=True, random_state=999)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**2. Train a non-regularized logistic regression model**\n",
    "\n",
    "Use parameters: penalty='none', solver='lbfgs', max_iter=200. Train on X_train, y_train."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "lg = ...\n",
    "# train model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<details>\n",
    "<summary>Solution</summary>\n",
    "<pre>\n",
    "lg = LogisticRegression(penalty='none', solver='lbfgs', max_iter=200)\n",
    "lg.fit(X_train, y_train)\n",
    "</pre>\n",
    "</details>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**3. Compare metrics for training and test set**\n",
    "\n",
    "See [model assessment](https://github.com/parrt/msds621/blob/master/lectures/model-assessment.pdf) for details on the log loss metric, but the score of zero means \"no loss\" (i.e., perfect score) and it is an unbounded positive loss from there, depending on how bad the model is. Log loss penalizes models that are confident in the wrong answer."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(f\"Accuracy: Train {100*lg.score(X_train, y_train):.0f}%, Test {100*lg.score(X_test, y_test):.0f}%\")\n",
    "y_proba_train = lg.predict_proba(X_train)[:,1]\n",
    "y_proba_test = lg.predict_proba(X_test)[:,1]\n",
    "print(f\"Log loss: Train {log_loss(y_train, y_proba_train):.2f}, Test {log_loss(y_test, y_proba_test):.2f}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q.**  What are some possible reasons for the difference in training and testing scores?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<details>\n",
    "<summary>Solution</summary>\n",
    "The model could be simply too weak. Our validation set could, by chance, look very different from the training set.  Outliers in the training or validation set can cause such differences.\n",
    "</details>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**4. Extract $\\beta_1, ..., \\beta_p$ and count how many are close to 0**\n",
    "\n",
    "Note: `sum(np.abs(x) < 1e-5)` is a decent way to check for values of `x` close to zero but not necessarily zero.  There is also `numpy.isclose()` but that is too strict (requires numbers to be really close to zero) for this exercise."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "lg_beta = ...\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<details>\n",
    "<summary>Solution</summary>\n",
    "<pre>\n",
    "lg_beta = lg.coef_[0]\n",
    "sum(np.abs(lg_beta) < 1e-5) # how many close to 0?\n",
    "</pre>\n",
    "</details>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**5. Plot the coefficient index versus the value**\n",
    "\n",
    "The plot should look something like:\n",
    "\n",
    "<img src=\"wine-ols.png\" width=\"200\">"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The following function is a handy way to plot the coefficients."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plotbeta(beta, which, yrange=(-20_000, 20_000),fontsize=10, xlabel=True, ylabel=True, tick_format='{x:.1f}', ax=None):\n",
    "    if ax is None:\n",
    "        fig,ax = plt.subplots(1,1,figsize=(4.5,2.5))\n",
    "    ax.bar(range(len(beta)),beta)\n",
    "    if xlabel:\n",
    "        ax.set_xlabel(\"Coefficient $\\\\beta_i$ for $i \\\\geq 1$\", fontsize=fontsize)\n",
    "    if ylabel:\n",
    "        ax.set_ylabel(\"Coefficient value\", fontsize=fontsize)\n",
    "    if yrange is not None:\n",
    "        ax.set_ylim(*yrange)\n",
    "    ax.tick_params(axis='both', which='major', labelsize=fontsize)\n",
    "    plt.gca().yaxis.set_major_formatter(mpl.ticker.StrMethodFormatter(tick_format))\n",
    "    ax.set_title(f\"{which} $\\\\overline{{\\\\beta}}$={np.mean(beta):.2f}, $\\\\sigma(\\\\beta)$={np.std(beta):.2f}\", fontsize=fontsize)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plotbeta(lg_beta, \"Non-regularized\", yrange=None)\n",
    "#plt.tight_layout(); plt.savefig(\"wine-ols.png\",dpi=150,bbox_inches=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**6. Normalize the variables, retrain, check how many coefficients are close to 0**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "normalize(X_train)\n",
    "normalize(X_test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "lg = LogisticRegression(penalty='none', solver='lbfgs')\n",
    "lg.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "lg_beta = lg.coef_[0]\n",
    "sum(np.abs(lg_beta) < 1e-5) # how many close to 0?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**7. Compare metrics for training and test set**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(f\"Accuracy: Train {100*lg.score(X_train, y_train):.0f}%, Test {100*lg.score(X_test, y_test):.0f}%\")\n",
    "y_proba_train = lg.predict_proba(X_train)[:,1]\n",
    "y_proba_test = lg.predict_proba(X_test)[:,1]\n",
    "print(f\"Log loss: Train {log_loss(y_train, y_proba_train):.2f}, Test {log_loss(y_test, y_proba_test):.2f}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q.**  Without the `max_iter` arg > 100 the logistic regression model for unnormalized data does not converge on a solution. Why does the normalized data converge faster?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<details>\n",
    "<summary>Solution</summary>\n",
    "The loss function contours after normalization are more spherical, which can improve the speed of convergence. It's also the case that zero-centering ensures that we have some positive and some negative gradients; otherwise all gradient descent motion has to move either all positively or all negatively, leading to zigzagging. More on this when we go into gradient descent.\n",
    "</details>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q.**  Why does the test set score improve just by normalizing variables?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<details>\n",
    "<summary>Solution</summary>\n",
    "This could have something to do with finite precision arithmetic or simply that it's easier to find the minimum for a more spherical loss function surface.\n",
    "</details>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**7. Plot the coefficient index versus the value again after normalizing variables**\n",
    "\n",
    "The plot should look something like:\n",
    "\n",
    "<img src=\"wine-ols-norm.png\" width=\"200\">"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plotbeta(lg_beta, \"Non-regularized normed\", yrange=None)\n",
    "#plt.tight_layout(); plt.savefig(\"wine-ols-norm.png\",dpi=150,bbox_inches=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q.** Why is the scale of the coefficients different after normalizing the variables?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<details>\n",
    "<summary>Solution</summary>\n",
    "The scale the coefficients is a function of the X variable ranges so we would expect normalized variables to  yield smaller coefficients, unless of course the range of the X variables was already small.\n",
    "</details>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q.** Why does the shape of the coefficient graph look different for the normalized data?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<details>\n",
    "<summary>Solution</summary>\n",
    "When some coefficients get very large, due to unnormalized data, coefficients with smaller ranges might get less emphasis. But, when all data is in the same range, we get a much more accurate picture of importance of the coefficients.\n",
    "</details>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## L1 Regularization"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\"\n",
    "sklearn says LogisticRegression arg C is \"Inverse of regularization strength...\n",
    "smaller values specify stronger regularization\"\n",
    "\"\"\"\n",
    "lmbda=.1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "lg = LogisticRegression(C=1/lmbda, penalty='l1', solver='liblinear', max_iter=1000)\n",
    "lg.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "lg_beta = lg.coef_[0]\n",
    "sum(np.abs(lg_beta) < 1e-5) # how many close to 0?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q.**  What does it mean for our model when 4 coefficients go to zero?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<details>\n",
    "<summary>Solution</summary>\n",
    "This gives us fewer variables in our model, which leads to a simpler and more efficient model. All else being equal, a model with fewer parameters often generalizes better.\n",
    "</details>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**2. Compare metrics for training and test set**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(f\"Accuracy: Train {100*lg.score(X_train, y_train):.0f}%, Test {100*lg.score(X_test, y_test):.0f}%\")\n",
    "y_proba_train = lg.predict_proba(X_train)[:,1]\n",
    "y_proba_test = lg.predict_proba(X_test)[:,1]\n",
    "print(f\"Log loss: Train {log_loss(y_train, y_proba_train):.2f}, Test {log_loss(y_test, y_proba_test):.2f}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**1. Plot the coefficient index versus the value**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plotbeta(lg_beta, \"Lasso\", yrange=None, tick_format='{x:.1f}')\n",
    "#plt.tight_layout(); plt.savefig(\"wine-l1.png\",dpi=150,bbox_inches=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q.** Compare the L1 regularized coefficients with the normalized OLS coefficients?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<details>\n",
    "<summary>Solution</summary>\n",
    "The shape of the coefficients appears to be the same but the scale is smaller for regularized coefficients. More importantly, 4 of the coefficients have gone to zero, which leads to a simpler and more efficient model.\n",
    "</details>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q.** Compare the accuracy of the regularized model to the normalized OLS model?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<details>\n",
    "<summary>Solution</summary>\n",
    "The classification error rate, the accuracy, is the same at 96% and the log loss is about the same. Here, the benefit is likely the fact that L1 kills some coefficients.</details>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Effect of $\\lambda$ on regularization and accuracy scores\n",
    "\n",
    "The goal of the following code snippets is to help you visualize how the $\\lambda$ regularization parameter affects model parameters and associated training and testing accuracy. There are a number of important questions to answer following the code snippets."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig,axes = plt.subplots(1,6,figsize=(11,1.8), sharey=True)\n",
    "for i,lmbda in enumerate([1e-5,.01,.1, 1, 10, 100]):\n",
    "    lg = LogisticRegression(C=1/lmbda, penalty='l1', solver='liblinear', max_iter=1000)\n",
    "    lg.fit(X_train, y_train)\n",
    "    accur = lg.score(X_train, y_train)\n",
    "    accurt = lg.score(X_test, y_test)\n",
    "    y_proba_train = lg.predict_proba(X_train)[:,1]\n",
    "    y_proba_test = lg.predict_proba(X_test)[:,1]\n",
    "    print(f\"lambda={lmbda:5}: Zeros {sum(np.abs(lg.coef_[0]) < 1e-5):3d}: Accuracy: Train {100*lg.score(X_train, y_train):3.0f}%, Test {100*lg.score(X_test, y_test):3.0f}%; Log loss: Train {log_loss(y_train, y_proba_train):.2f}, Test {log_loss(y_test, y_proba_test):.2f}\")\n",
    "    plotbeta(lg.coef_[0], f\"$\\\\lambda={lmbda}$\\n\", ax=axes[i], fontsize=9, xlabel=False, ylabel=i==0, yrange=None)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q.** Describe what is happening to the L1 coefficients as we increase lambda?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<details>\n",
    "<summary>Solution</summary>\n",
    "The average magnitude of the coefficients is changing and they are becoming a tighter group; the standard deviation is shrinking significantly. Many more coefficients are going to zero.\n",
    "</details>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q.** Why does the training accuracy and log loss go down as we increase lambda?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<details>\n",
    "<summary>Solution</summary>\n",
    "We are introducing bias, in exchange for increased generality, but clearly you can increase $\\lambda$ too much and kill accuracy even for the training set.\n",
    "</details>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q.** What $\\lambda$ value would you choose for regularization?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<details>\n",
    "<summary>Solution</summary>\n",
    "Based upon the tests above, $\\lambda=0.1$ seems a good choice because it gets the lowest log loss on the test set.\n",
    "</details>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## L2 Regularization"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "lmbda=0.01"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "lg = LogisticRegression(C=1/lmbda, penalty='l2', solver='liblinear', max_iter=1000)\n",
    "lg.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "lg_beta = lg.coef_[0]\n",
    "sum(np.abs(lg_beta) < 1e-5) # how many close to 0?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**2. Compare the R^2 on training and test set**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(f\"Accuracy: Train {100*lg.score(X_train, y_train):.0f}%, Test {100*lg.score(X_test, y_test):.0f}%\")\n",
    "y_proba_train = lg.predict_proba(X_train)[:,1]\n",
    "y_proba_test = lg.predict_proba(X_test)[:,1]\n",
    "print(f\"Log loss: Train {log_loss(y_train, y_proba_train):.2f}, Test {log_loss(y_test, y_proba_test):.2f}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**1. Plot the coefficient index versus the value**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plotbeta(lg_beta, \"Ridge\", yrange=None, tick_format='{x:.1f}')\n",
    "#plt.tight_layout(); plt.savefig(\"wine-l1.png\",dpi=150,bbox_inches=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Effect of $\\lambda$ on regularization and accuracy scores\n",
    "\n",
    "The goal of the following code snippets is to help you visualize how the $\\lambda$ regularization parameter affects model parameters and associated training and testing accuracy. There are a number of important questions to answer following the code snippets."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig,axes = plt.subplots(1,6,figsize=(11,1.8), sharey=True)\n",
    "for i,lmbda in enumerate([1e-5,.01,.1, 1, 10, 100]):\n",
    "    lg = LogisticRegression(C=1/lmbda, penalty='l2', solver='liblinear', max_iter=1000)\n",
    "    lg.fit(X_train, y_train)\n",
    "    accur = lg.score(X_train, y_train)\n",
    "    accurt = lg.score(X_test, y_test)\n",
    "    y_proba_train = lg.predict_proba(X_train)[:,1]\n",
    "    y_proba_test = lg.predict_proba(X_test)[:,1]\n",
    "    print(f\"lambda={lmbda:5}: Zeros {sum(np.abs(lg.coef_[0]) < 1e-5):3d}: Accuracy: Train {100*lg.score(X_train, y_train):3.0f}%, Test {100*lg.score(X_test, y_test):3.0f}%; Log loss: Train {log_loss(y_train, y_proba_train):.2f}, Test {log_loss(y_test, y_proba_test):.2f}\")\n",
    "    plotbeta(lg.coef_[0], f\"$\\\\lambda={lmbda}$\\n\", ax=axes[i], fontsize=9, xlabel=False, ylabel=i==0, yrange=None)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q.** Describe what is happening to the L2 coefficients as we increase lambda?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<details>\n",
    "<summary>Solution</summary>\n",
    "The average magnitude of the coefficients goes down and they become a tighter group. We don't get any  coefficients going to zero.\n",
    "</details>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q.** Characterize how the magnitudes of L1 and L2 coefficients differ."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<details>\n",
    "<summary>Solution</summary>\n",
    "   The L2 coefficients are in general much tighter group than the L1, even as we increase $\\lambda$. The L1 regularization drops many of the coefficients to zero for the same value of $\\lambda$.\n",
    "</details>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q.**  Which regularization method would you choose for this data set?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<details>\n",
    "<summary>Solution</summary>\n",
    "L2 regularization is a clear winner here. Not only is the accuracy on the test set 100%, but the log loss on the test set is smaller than any other technique.\n",
    "</details>"
   ]
  }
 ],
 "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.8.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}