15_ml.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Machine Learning in Python\n",
    "\n",
    "- Tools: scikit-learn (`sklearn`)\n",
    "    - Data Partitioning\n",
    "    - Feature selection\n",
    "    - Modeling: SVM\n",
    "    - Model Assessment\n",
    "    \n",
    " [![Open In Colab](https://colab.research.google.com/assets/colab-badge.svg)](https://colab.research.google.com/github/COGS108/Lectures-Fa24/blob/main/15_ml.ipynb)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For more reading on scikit-learn (`sklearn`) and machine learning in Python: https://scikit-learn.org/stable/index.html"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Machine Learning: General Steps\n",
    "\n",
    "1. Data Partitioning\n",
    "2. Feature Selection\n",
    "3. Model\n",
    "4. Model Assessment"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Setup"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# import ds/plotting packages\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "\n",
    "%matplotlib inline\n",
    "\n",
    "# Import nltk package \n",
    "import nltk\n",
    "from nltk.tokenize import word_tokenize\n",
    "from nltk.corpus import stopwords\n",
    "\n",
    "# import random for randomizing\n",
    "import random\n",
    "\n",
    "# ignore warnings\n",
    "import warnings\n",
    "warnings.filterwarnings('ignore')\n",
    "\n",
    "# scikit-learn imports\n",
    "from sklearn.svm import SVC\n",
    "from sklearn.feature_extraction.text import TfidfVectorizer\n",
    "from sklearn.metrics import classification_report, precision_recall_fscore_support, confusion_matrix, ConfusionMatrixDisplay"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Uncomment if you need to download the NLTK English tokenizer and the stopwords of all languages\n",
    "nltk.download('punkt')\n",
    "nltk.download('stopwords')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Example: Train-test split with penguin data\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# lets get some 🐧 data to work with\n",
    "penguins = sns.load_dataset(\"penguins\").dropna()\n",
    "\n",
    "sns.pairplot(penguins, hue='species')\n",
    "\n",
    "penguins"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.model_selection import train_test_split\n",
    "\n",
    "# here's our raw data\n",
    "X = penguins[['bill_length_mm', 'bill_depth_mm', \n",
    "                  'flipper_length_mm', 'body_mass_g']] #, 'sex']]\n",
    "y = penguins['species']\n",
    "\n",
    "# make a train test split with params\n",
    "# test_size=0.33, random_state=42\n",
    "# see sklearn docs!\n",
    "\n",
    "X_train, X_test, y_train, y_test = train_test_split(\n",
    "    X, y, test_size=0.33, random_state=42)\n",
    "\n",
    "print(X_train.shape, X_test.shape, y_train.shape, y_test.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# models we will use\n",
    "from sklearn.tree import DecisionTreeClassifier\n",
    "\n",
    "# here's our model!  look at the docs for options!!\n",
    "# right now we are going to use the defaults.\n",
    "model = DecisionTreeClassifier()\n",
    "\n",
    "model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model.fit(X_train, y_train )\n",
    "\n",
    "\n",
    "# how good will we do on training set error and test set error?\n",
    "yhat_train = model.predict(X_train)\n",
    "yhat_test = model.predict(X_test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.metrics import accuracy_score, classification_report\n",
    "from sklearn.metrics import ConfusionMatrixDisplay\n",
    "\n",
    "acc_train = accuracy_score(y_train, yhat_train)\n",
    "\n",
    "# calculate the test set error\n",
    "acc_test = accuracy_score(y_test, yhat_test)\n",
    "\n",
    "print(f'training set accuracy (n={y_train.shape[0]}): {acc_train:4.3f}')\n",
    "print(f'test set accuracy (n={y_test.shape[0]}): {acc_test:4.3f}')\n",
    "print()\n",
    "print('classificiation report on test set performance')\n",
    "print(classification_report(y_test, yhat_test))\n",
    "print()\n",
    "\n",
    "with sns.axes_style('white'):\n",
    "    ConfusionMatrixDisplay.from_predictions(y_test, yhat_test)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# What if I want to do model selection using a validation set?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Defines ratios, w.r.t. whole dataset.\n",
    "ratio_train = 0.8\n",
    "ratio_val = 0.1\n",
    "ratio_test = 0.1\n",
    "\n",
    "\n",
    "# here's our raw data\n",
    "X = penguins[['bill_length_mm', 'bill_depth_mm', \n",
    "                  'flipper_length_mm', 'body_mass_g']] #, 'sex']]\n",
    "y = penguins['species']\n",
    "\n",
    "# Produces test split.\n",
    "x_remaining, x_test, y_remaining, y_test = train_test_split(\n",
    "    X, y, test_size=ratio_test)\n",
    "\n",
    "# Adjusts val ratio, w.r.t. remaining dataset.\n",
    "ratio_remaining = 1 - ratio_test\n",
    "ratio_val_adjusted = ratio_val / ratio_remaining\n",
    "\n",
    "# Produces train and val splits.\n",
    "x_train, x_val, y_train, y_val = train_test_split(\n",
    "    x_remaining, y_remaining, test_size=ratio_val_adjusted)\n",
    "\n",
    "print(x_train.shape, x_val.shape, x_test.shape, y_train.shape, y_val.shape, y_test.shape)\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## now you can...\n",
    "- setup 3 decision tree models with different values of `max_depth` in the range `1, 2, 3`\n",
    "- fit all 3 models on `x_train, y_train`\n",
    "- get `accuracy_score(x_val, y_val)`\n",
    "- pick the best of the models, see its `accuracy_score(x_test, y_test)`\n",
    "\n",
    "exercise left to you to figure out for yourself!\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Example: transforming (pre-processing) data for better performance\n",
    "One thing I was wondering, is it possible to increase the accuracy of the predictions by including sex.  A bit of domain knowledge is helpful here... many species have different sized males and females, and you could easily see that impacting prediction accuracy, where males of smaller species might have measurements similar to females of larger species, causing misclassification."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# make a combo category\n",
    "penguins['species_x_sex'] = penguins['species']+', '+penguins['sex']\n",
    " \n",
    "# this makes color sequence light/dark of the same color\n",
    "# to encode sex for a given species, see seaborn docs\n",
    "with sns.color_palette(\"Paired\"):\n",
    "\n",
    "    # pariplot shows all pairs of real valued variables\n",
    "    sns.pairplot(penguins, hue='species_x_sex', \n",
    "                 hue_order=[ # forces the right order of colors \n",
    "                    'Adelie, Male', 'Adelie, Female',\n",
    "                    'Chinstrap, Male', 'Chinstrap, Female',\n",
    "                    'Gentoo, Male', 'Gentoo, Female']\n",
    "                );"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Yep, it sure looks like there is indeed a systematic sex based difference in size across all 3 species.  so lets do it again, this time we will try to improve our classification accuracy by including sex as a feature.\n",
    "\n",
    "Now one thing... sex is a categorical feature, and so it needs transformation. And while its a binary feature in this dataset, we will just go ahead and one-hot encode it for safety.  In fact, while some sklearn algorithms will handle categorical inputs with ease, decision tree classifier is one of the ones that requires us to one-hot a categorical feature before it will work right.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# here's our raw data, this time include sex!\n",
    "X = penguins[['bill_length_mm', 'bill_depth_mm', \n",
    "                  'flipper_length_mm', 'body_mass_g', 'sex']]\n",
    "y = penguins['species']\n",
    "\n",
    "# here's our training set and test set\n",
    "X_train, X_test, y_train, y_test = train_test_split(\n",
    "    X, y, test_size=0.33, random_state=42)\n",
    "\n",
    "print(X_train.shape, X_test.shape, y_train.shape, y_test.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.preprocessing import OneHotEncoder, StandardScaler\n",
    "\n",
    "# transformers\n",
    "onehot = OneHotEncoder() # use me for categorical data\n",
    "scaler = StandardScaler() # for some algos, you need me for numerical data\n",
    "                          # for other algos, you don't need it!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Whenever we do any kind of preprocessing \n",
    "# (e.g. transformation or normalization)\n",
    "# to our data it is CRITICAL to fit the preprocessor on training data\n",
    "# NOT all the data.  If you fit all the data then you are leaking\n",
    "# information from the test set into training set, and therefore\n",
    "# you are allowing overfitting\n",
    "\n",
    "onehot.fit( X_train[['sex']] )\n",
    "\n",
    "# the to_frame() thing is only necessary becuase onehot expects a\n",
    "# 2d array, and a single column is a 1d array.  if we have multiple\n",
    "# columns to onehot we wouldn't need to_frame() and it wouldn't \n",
    "# exist anyway as it would already be a 2d dataframe"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# add new one hot encoded gender categories to datasets\n",
    "# ends up as a numpy sparse matrix, so we need to densify it\n",
    "# to assign it back into the dataframe\n",
    "X_train.loc[:, onehot.categories_[0]] = onehot.transform(X_train[['sex']]).todense()\n",
    "X_test.loc[:, onehot.categories_[0]] = onehot.transform(X_test[['sex']]).todense()\n",
    "\n",
    "# drop the old text label category \"sex\"\n",
    "X_train.drop('sex', axis='columns', inplace=True)\n",
    "X_test.drop('sex', axis='columns', inplace=True)\n",
    "\n",
    "X_train"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model.fit(X_train,y_train);\n",
    "\n",
    "yhat_train = model.predict(X_train)\n",
    "yhat_test = model.predict(X_test)\n",
    "\n",
    "acc_train = accuracy_score(y_train, yhat_train)\n",
    "acc_test = accuracy_score(y_test, yhat_test)\n",
    "\n",
    "print(f'training set accuracy (n={y_train.shape[0]}): {acc_train:4.3f}')\n",
    "print(f'test set accuracy (n={y_test.shape[0]}): {acc_test:4.3f}')\n",
    "print()\n",
    "print('classificiation report on test set performance')\n",
    "print(classification_report(y_test, yhat_test))\n",
    "print()\n",
    "\n",
    "# use this `with` to avoid ugly gridlines in the confusion matrix plot\n",
    "# that an earlier seaborn style statement would have put there\n",
    "# this with is not necessary if you havaen't been using seaborn above\n",
    "with sns.axes_style('white'):\n",
    "    ConfusionMatrixDisplay.from_predictions(y_test, yhat_test)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Example: Cross validation instead of train/test split\n",
    "\n",
    "What if we want to use cross-validation to do evaluation generalization?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.model_selection import cross_val_score\n",
    "from sklearn.svm import LinearSVC \n",
    "\n",
    "y_p = penguins['species'] # predict species from the measurements\n",
    "X_p = penguins.drop(columns=['sex','species','island','species_x_sex']) # get rid of categorical variables and the label\n",
    "   \n",
    "    \n",
    "results = cross_val_score( LinearSVC(), # here we are using default values for SVC!!! you need to look at the docs!!! \n",
    "                           X_p, y_p, cv=5)\n",
    "\n",
    "print(results)\n",
    "print( np.mean(results))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Example: Model selection and evaluation using cross-validation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.neighbors import KNeighborsClassifier\n",
    "\n",
    "from sklearn.model_selection import train_test_split\n",
    "\n",
    "# here's our raw data\n",
    "X = penguins[['bill_length_mm', 'bill_depth_mm', \n",
    "                  'flipper_length_mm', 'body_mass_g']]\n",
    "y = penguins['species']\n",
    "\n",
    "# make a train test split with params\n",
    "# test_size=0.33, random_state=42\n",
    "# see sklearn docs!\n",
    "\n",
    "X_train, X_test, y_train, y_test = train_test_split(\n",
    "    X, y, test_size=0.33, random_state=42)\n",
    "\n",
    "\n",
    "\n",
    "# how many neighbors is best?\n",
    "for K in [1, 3, 5, 7, 9, 11, 13, 15, 17, 19]:\n",
    "    knn = KNeighborsClassifier(n_neighbors=K)\n",
    "    results = cross_val_score( knn, X_train, y_train, cv=5)\n",
    "    print('For ',K,' neighbors validation accuracy is ',np.mean(results)) \n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "best = KNeighborsClassifier(n_neighbors=1)\n",
    "best.fit(X_train, y_train)\n",
    "\n",
    "yhat_train = best.predict(X_train)\n",
    "yhat_test = best.predict(X_test)\n",
    "\n",
    "acc_train = accuracy_score(y_train, yhat_train)\n",
    "acc_test = accuracy_score(y_test, yhat_test)\n",
    "\n",
    "print(f'training set accuracy (n={y_train.shape[0]}): {acc_train:4.3f}')\n",
    "print(f'test set accuracy (n={y_test.shape[0]}): {acc_test:4.3f}')\n",
    "print()\n",
    "print('classificiation report on test set performance')\n",
    "print(classification_report(y_test, yhat_test))\n",
    "print()\n",
    "\n",
    "# use this `with` to avoid ugly gridlines in the confusion matrix plot\n",
    "# that an earlier seaborn style statement would have put there\n",
    "# this with is not necessary if you havaen't been using seaborn above\n",
    "with sns.axes_style('white'):\n",
    "    ConfusionMatrixDisplay.from_predictions(y_test, yhat_test)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Example: using TF-IDF from class comments to predict whether the comment is the most/lest liked \n",
    "\n",
    "## Data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Student responses on COGS 108 Mid-course survey to the following two questions: \n",
    "\n",
    "- What have you enjoyed MOST about COGS 108 so far? Please explain.\n",
    "- What have you enjoyed LEAST about COGS 108 so far? Please explain."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# read data in\n",
    "# 1 = most; 0 = least\n",
    "df = pd.read_csv('https://raw.githubusercontent.com/shanellis/datasets/master/COGS108_ml.csv', encoding=\"ISO-8859-1\")\n",
    "df.tail()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# randomly sort data frame\n",
    "df = df.sample(frac=1, random_state=200).reset_index(drop=True)\n",
    "df.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Randomly sorted data frame:\n",
    "- for selection of training and test set\n",
    "- will be approximately balanced between outcomes in each"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# see how much data we're working with\n",
    "df.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df.groupby('quarter')['response'].count()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Train,  Validate, Test\n",
    "\n",
    "- Training set is what we use for creating the model\n",
    "- Validation set is what we use for fine tuning the model -- checking which hyper-parameters are best\n",
    "- Test set is what we use to figure out how good the model will do next year with new data arriving\n",
    "\n",
    "#### Method (1) - Manual split, split by quarter\n",
    "\n",
    "- We'll test the model on one quarter's responses\n",
    "- We'll validate the model on another quarter's responses\n",
    "- We'll train the model on everything else\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "\n",
    "## Validation\n",
    "df_validation = df[df['quarter']=='sp20']\n",
    "              \n",
    "## Test\n",
    "df_test = df[df['quarter']=='sp19']\n",
    "\n",
    "## Train\n",
    "df_train = df[~df['quarter'].isin(['sp19','sp20'])]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(df_train.shape, df_validation.shape, df_test.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Clicker Question #1\n",
    "\n",
    "How well do you think we will be able to predict whether a student comment is a response to what they liked most vs what they liked least in our ***training*** dataset?\n",
    "\n",
    "- A) Accuracy ~0%\n",
    "- B) Accuracy ~25%\n",
    "- C) Accuracy ~50%\n",
    "- D) Accuracy ~75%\n",
    "- E) Accuracy ~100%"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Clicker Question #2\n",
    "\n",
    "How well do you think we will be able to predict whether a student comment is a response to what they liked most vs what they liked least in our ***validation*** dataset?\n",
    "\n",
    "- A) Accuracy ~0%\n",
    "- B) Accuracy ~25%\n",
    "- C) Accuracy ~50%\n",
    "- D) Accuracy ~75%\n",
    "- E) Accuracy ~100%"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Thinking this through\n",
    "\n",
    "Is this a good idea to split up by quarter?  What could happen here?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Method (2) - sklearn split, randomize the data\n",
    "- Training set 60% of the data randomly mixed across quarters\n",
    "- Validation set 20% of the data randomly mixed across quarters\n",
    "- Test set 20% of the data randomly mixed across quarters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.model_selection import train_test_split\n",
    "\n",
    "X_trainvalid, X_test, Y_trainvalid, Y_test = train_test_split(\n",
    "    df['response'], df['most_least'], test_size=0.2, random_state=42)\n",
    "\n",
    "X_train, X_valid, Y_train, Y_valid = train_test_split(\n",
    "    X_trainvalid, Y_trainvalid, test_size=0.25\n",
    "\n",
    ")\n",
    "\n",
    "print(X_train.shape, X_valid.shape, X_test.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Prediction Task: \n",
    "\n",
    "**Classify text from students as 'most liked' or 'least liked'**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 11 Steps to Prediction:\n",
    "\n",
    "1. Specify parameters for TF-IDF calculation\n",
    "2. Calculate TF-IDF from text input (predictors)\n",
    "3. Extract most or least (outcome)\n",
    "4. Specify how data will be partitioned\n",
    "5. Partition the data\n",
    "6. Train model\n",
    "7. Predict in training\n",
    "8. Predict in testing\n",
    "9. Assess accuracy in training\n",
    "10. Assess accuracy in test set\n",
    "11. IF NEEDED: Model selection + test process that replaces steps 6-10 "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Data Processing"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 1: Determine how you'll convert a collection of raw documents to a matrix of TF-IDF features."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create vectorizer & specify parameters\n",
    "tfidf = TfidfVectorizer(sublinear_tf=True, #apply sublinear TF scaling\n",
    "                        analyzer='word',   #specify tokenizer\n",
    "                        max_features=500, # specify max # of features to include\n",
    "                        tokenizer=word_tokenize)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* sublinear TF scaling - replaces term frequency (TF) with $1 + log(TF)$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 2: Generate matrix of TF-IDF features."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Learn vocabulary and idf, return term-document matrix.\n",
    "# return an array;our predictor\n",
    "tfidf_X = tfidf.fit_transform(X_train).toarray()\n",
    "\n",
    "# take a look at the output\n",
    "print(tfidf_X.shape)\n",
    "\n",
    "print(\"min: \" , np.min(tfidf_X), '\\n',\n",
    "      \"mean: \", np.mean(tfidf_X), '\\n',\n",
    "      \"max: \",  np.max(tfidf_X))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "tfidf_X"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "## get IDF to visualize\n",
    "idf = tfidf.idf_\n",
    "rr = dict(zip(tfidf.fit(X_train).get_feature_names(), idf))\n",
    "\n",
    "token_weight = pd.DataFrame.from_dict(rr, orient='index').reset_index()\n",
    "token_weight.columns=('token','weight')\n",
    "token_weight = token_weight.sort_values(by='weight', ascending=False)\n",
    "token_weight"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sns.barplot(x='token', \n",
    "            y='weight', \n",
    "            data=token_weight[0:10], \n",
    "            color=\"gray\")            \n",
    "plt.title(\"Inverse Document Frequency(IDF) per token\")\n",
    "fig = plt.gcf()\n",
    "fig.set_size_inches(15,5);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 3: Extract outcome variable"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# specify outcome variable\n",
    "tfidf_Y = np.array(Y_train)\n",
    "tfidf_Y"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "tfidf_X"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "tfidf_X.shape, tfidf_Y.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Clicker Question #3\n",
    "\n",
    "Looking at the code above and thinking about what we've done so far in this analysis, what is stored in `tfidf_Y`?\n",
    "\n",
    "- A) predictor variable - training data\n",
    "- B) outcome variable - training data\n",
    "- C) predictor variable - test data\n",
    "- D) outcome variable - test data\n",
    "- E) validation DataFrame"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Avoiding data leakage with feature selection\n",
    "\n",
    "NEVER do feature selection or feature transformation fitting on ALL of the data.\n",
    "\n",
    "You don't want information about your validation and test set inside your training set.\n",
    "\n",
    "If you fit your features on the whole dataset you have implictly told the training set information about the vocuabulary and frequency of things that aren't in the training set! BAD!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "tfidf_train = tfidf.fit_transform(X_train)\n",
    "tfidf_valid = tfidf.transform(X_valid)\n",
    "tfidf_test = tfidf.transform(X_test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "tfidf_train.shape, tfidf_valid.shape, tfidf_test.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### SVM: Support Vector Machines\n",
    "\n",
    "- simple & interpretable machine learning model\n",
    "- based in linear regression\n",
    "- classification task\n",
    "- supervized\n",
    "    - input: labeled training data\n",
    "    - model determines hyperplane that best discriminates between categories"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### SVM: Tuning Parameters\n",
    "- **regularization** parameter\n",
    "    - can determine how this line is drawn\n",
    "    - can increase accuracy of prediction\n",
    "    - can lead to overfitting of the data\n",
    "- **kernel** parameter\n",
    "    - specifies how to model & transform data\n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For more reading on SVMs using `sklearn`: https://scikit-learn.org/stable/modules/svm.html\n",
    "\n",
    "OR just scroll to the very end of this notebook where we get a bit crazy with them "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Model Generation\n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 6: Generate and train the model."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# uncomment to read documentation for model\n",
    "SVC?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# function we'll use to run the model\n",
    "\n",
    "svm_model = SVC(kernel='linear')\n",
    "type(svm_model)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# train model\n",
    "svm_model = svm_model.fit(tfidf_train, Y_train)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Training Data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 7: Predict in the training data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# predict on training\n",
    "df_predicted_train_Y = svm_model.predict(tfidf_train)\n",
    "\n",
    "print(df_predicted_train_Y[0:5])\n",
    "print(Y_train[0:5])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# see how many were predicted most vs. least\n",
    "pd.Series(df_predicted_train_Y).value_counts()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Testing Data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 8: Predict in the testing data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# predict on training\n",
    "df_predicted_test_Y = svm_model.predict(tfidf_test)\n",
    "print(df_predicted_test_Y[0:5])\n",
    "print(Y_test[0:5])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# see how many were predicted most vs. least\n",
    "pd.Series(df_predicted_test_Y).value_counts()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Accuracy Assessment\n",
    "\n",
    "- RMSE (continuous)\n",
    "- Accuracy, Sensitivity, Specificity, AUC\n",
    "    - TP, TN, FP, FN"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![confusion matrix](https://miro.medium.com/max/924/1*7EYylA6XlXSGBCF77j_rOA.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![sensitivity recell](https://miro.medium.com/max/878/1*Ub0nZTXYT8MxLzrz0P7jPA.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Accuracy** - What % were predicted correctly?  \n",
    "**Sensitivity (Recall)** - Of those that were positives, what % were predicted to be positive?  ; $\\frac {TP}{(TP + FN)}$  \n",
    "**Specificity** - Of those that were actually negatives, what % were predicted to be negative?  $\\frac {TN}{(TN + FP)}$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Precision (Positive Predictive Value, PPV)** = $\\frac {TP}{(TP + FP)}$\n",
    "\n",
    "- probability that predicted positive truly is positive"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Training Data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 9: Assess accuracy in training data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(classification_report(Y_train, df_predicted_train_Y))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# where 'support' comes from\n",
    "pd.Series(Y_train).value_counts()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cm = confusion_matrix(Y_train, df_predicted_train_Y)\n",
    "disp = ConfusionMatrixDisplay(cm)\n",
    "disp.plot();\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**support** - the number of occurrences of each class  \n",
    "**precision (PPV)** - ability of the classifier not to label a positive sample as negative  \n",
    "**recall (sensitivity)** - ability of the classifer to find all the positive samples\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**f1-score** - weighted harmonic mean of the precision and recall; score reaches its best value at 1 and worst score at 0  \n",
    "**macro average** - averaging the unweighted mean per label  \n",
    "**weighted average** - averaging the support-weighted mean per label  \n",
    "**micro average** - averaging the total true positives, false negatives and false positives"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Testing Data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 10 or 11: Assess accuracy in testing data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(classification_report(Y_test, df_predicted_test_Y))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cm = confusion_matrix(Y_test, df_predicted_test_Y)\n",
    "disp = ConfusionMatrixDisplay(cm)\n",
    "disp.plot();"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Clicker Question #4\n",
    "\n",
    "Given this output, would you use this model to predict whether or not text was something someone liked or disliked about COGS 108?\n",
    "\n",
    "- A) Yes\n",
    "- B) No\n",
    "- C) Unsure"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Doing model selection with the validation data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Actually, usually BEFORE you test, you'll want to play with hyperparameters in validation to get the best result.  Fine tuning your model.\n",
    "\n",
    "\n",
    "So consider this an optional step, that usually comes before test if your model has a set of knobs and dials to twiddle.  \n",
    "\n",
    "We only did test first above to build some intuition for you."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# play with the parameter C until it looks good\n",
    "svm_model = SVC(kernel='linear', C=10)\n",
    "svm_model.fit(tfidf_train, Y_train)\n",
    "df_predicted_validation_Y = svm_model.predict(tfidf_valid)\n",
    "\n",
    "# assess accuracy\n",
    "cm = confusion_matrix(Y_valid, df_predicted_validation_Y)\n",
    "disp = ConfusionMatrixDisplay(cm)\n",
    "disp.plot();\n",
    "\n",
    "print(classification_report(Y_valid, df_predicted_validation_Y))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#assuming we left the cell above with our BEST model setting\n",
    "# lets refit it with all of train+valid sets\n",
    "svm_model.fit(np.concatenate([tfidf_train.todense(), tfidf_valid.todense()]), \n",
    "              np.concatenate([Y_train, Y_valid]))\n",
    "\n",
    "# get test set predictions\n",
    "df_predicted_test_Y = svm_model.predict(tfidf_test.todense())\n",
    "\n",
    "print(classification_report(Y_test, df_predicted_test_Y))\n",
    "\n",
    "cm = confusion_matrix(Y_test, df_predicted_test_Y)\n",
    "disp = ConfusionMatrixDisplay(cm)\n",
    "disp.plot();"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Summary"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "1. 80:20 Partition\n",
    "2. Specified TF-IDF as predictor and most/least (0,1) as outcome\n",
    "3. Trained SVM linear classifier\n",
    "4. Built model on Training data\n",
    "5. Predicted in training data and on testing data\n",
    "6. Assessed overall accuracy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Approaches For Improvement?\n",
    "\n",
    "- Data Cleaning/Stemming\n",
    "- Different Tuning Parameters?\n",
    "- Cross-Validation?\n",
    "- Train/Test on all data OR Train only Sp20/Fa20/Wi21 (remote quarters)\n",
    "- Different Model?\n",
    "\n",
    "To start with use sklearn the EASY way... don't try to hand code things like Train-Test split or K-Fold cross-validation... read the docs and learn how to use these wonderful functions https://scikit-learn.org/stable/modules/classes.html#module-sklearn.model_selection"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.model_selection import cross_val_score\n",
    "results = cross_val_score(svm_model, \n",
    "                         X=np.concatenate([tfidf_train.todense(), tfidf_valid.todense()]), \n",
    "                         y=np.concatenate([Y_train, Y_valid]),\n",
    "                         cv=5) \n",
    "import numpy as np\n",
    "np.mean(results)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Feature engineering\n",
    "One thing we didn't talk about is how you need to engineer features for NON-text kinds of inputs\n",
    "\n",
    "In general this is a huge topic that I will try condense for you into two bullet points.\n",
    "\n",
    "- If your data has different kinds of continuous variables, and especially if some of those variables tend to be big numbers (distance to the moon in meters) and some tend to be small numbers (the earths mass as a fraction of the earth's mass), then you need to STANDARDIZE those two variables to have roughly the same magnitude. If you don't then the ML system may think that the large magnitude variable is more predictive than the small magnitude variable.  See this link https://scikit-learn.org/stable/modules/preprocessing.html#preprocessing-scaler\n",
    "- If your data has categorical variables that do not have a natural order to them (red is not greater than blue), then you need to ONE-HOT encode those categorical variables.  Otherwise your ML system will think that red is greater than blue and that can make for some bad predictions.  See this link https://scikit-learn.org/stable/modules/preprocessing.html#preprocessing-categorical-features\n",
    "\n",
    "OK lets talk about how to do things with the non text input..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import seaborn as sns\n",
    "import pandas as pd \n",
    "\n",
    "%matplotlib inline\n",
    "\n",
    "# lets get some 🐧 data to work with\n",
    "df = sns.load_dataset('penguins').dropna()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sns.pairplot(df, hue='species');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Hmm... lets take a close look at body mass vs bill length.... that seems really promising!!\n",
    "\n",
    "Let's imagine what it would be like to use k-Nearest Neighbors on only those two inputs... it might be enough to make things predict well!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# what we want to predict\n",
    "y = df['species']\n",
    "# the data we will predict from\n",
    "X = df[['bill_length_mm','body_mass_g']]\n",
    "\n",
    "X"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "UGH!!! the vector difference $\\mathbf{x}_{new} -  \\mathbf{x}_{i}$ is going to be absolutely dominated by body mass... it's two orders of magnitude bigger than bill length!\n",
    "\n",
    "We may not want that!  From the scatterplot it looks like both variables are equally contibuting to understanding species.\n",
    "\n",
    "We need to normalize the data so each variable is a z-scored version of itself... that means both variables will now be on roughly the same order of magnitued (about -3 to 3).  There are in fact many diffrerent kinds of normalization that could work here!  Take a look at this list: https://scikit-learn.org/stable/modules/classes.html#module-sklearn.preprocessing"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.preprocessing import StandardScaler\n",
    "\n",
    "preproc = StandardScaler()\n",
    "\n",
    "preproc.fit(X)\n",
    "scaled_data = preproc.transform(X)\n",
    "X['bill_length_mm'] = scaled_data[:,0]\n",
    "X['body_mass_g'] = scaled_data[:,1]\n",
    "\n",
    "X\n",
    "\n",
    "# some things to note here... we could just use scaled_data instead of stuffing\n",
    "# the results back into the dataframe X... I'm just doing this so we can keep\n",
    "# the column labels intact... scaled_data is a numpy ndarray without those labels\n",
    "# sklearn will accept numpy arrays or pandas dataframes/series or any other iterable \n",
    "# as input, so use whatever you want"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The thing we decided to learn above is strictly numerical.  And preprocessing the way we just described is only good for numerical data.\n",
    "\n",
    "What if we had been interested in categorical data? For categorical data inputs we usually wish to use a one-hot encoding.\n",
    "\n",
    "One-hotting a single categorical variable that has N possible categories makes a vector representation that is N long. Because one-hot assumes each category is mutually exclusive with every other category, only one of those N category holders is allowed to be non-zero.\n",
    "\n",
    "Why use one-hot encoding?  This eliminates unwanted \"fake relationships\" in categorical data that a different encoding might cause an ML system to learn.\n",
    "\n",
    "Imagine we have categories like (red, blue, green). There is no ordinal relationship here... red is not > blue, nor is blue > green.  But if we turned (red, blue, green) into a representation like (3, 2, 1) that is exactly what our ML system would learn. Even worse, it would be more likely to confuse blue and green than red and green... because blue-green is smaller than red-green.\n",
    "\n",
    "One hot eliminates this because the 3 long vectors representing each color [1, 0, 0], [0, 1,0], and [0, 0, 0] are orthogonal (unrelated)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.preprocessing import OneHotEncoder\n",
    "\n",
    "preproc = OneHotEncoder()\n",
    "\n",
    "onehot = preproc.fit_transform(df[['island','sex']])\n",
    "print(df[['island','sex']])\n",
    "print(onehot.todense())\n",
    "\n",
    "# there are 3 different islands, so 3 columns get used to encode which island\n",
    "# only 1 of those 3 columns can be \"hot\" (non-zero) for each penguin\n",
    "# likewise there are 2 different sexes recorded here, so 2 columns get used"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "anyway we're not going to use categorical data right now, so lets stick to the numerical data we did above"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.neighbors import KNeighborsClassifier\n",
    "from sklearn.metrics import classification_report, ConfusionMatrixDisplay, confusion_matrix\n",
    "\n",
    "pclassifier = KNeighborsClassifier(n_neighbors=15)\n",
    "\n",
    "pclassifier.fit(X,y)\n",
    "predicted = pclassifier.predict(X)\n",
    "\n",
    "print(classification_report(y,predicted))\n",
    "\n",
    "cm = confusion_matrix(y,predicted)\n",
    "disp = ConfusionMatrixDisplay(cm)\n",
    "disp.plot();"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This result shows us that with just 2 variables we can get to 95% accuracy!\n",
    "\n",
    "BUT note this is TRAINING SET accuracy... something of limited use.\n",
    "\n",
    "Let's figure out how to make a training set/test set\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.model_selection import train_test_split\n",
    "\n",
    "X_train, X_test, y_train, y_test = train_test_split( X, y, test_size=0.3, random_state=777)\n",
    "\n",
    "print('training/test set sizes')\n",
    "print(X_train.shape, X_test.shape, y_train.shape, y_test.shape)\n",
    "print()\n",
    "\n",
    "pclassifier = KNeighborsClassifier(n_neighbors=15)\n",
    "\n",
    "pclassifier.fit(X_train, y_train)\n",
    "predicted_train = pclassifier.predict(X_train)\n",
    "predicted_test = pclassifier.predict(X_test)\n",
    "\n",
    "print('training set performance')\n",
    "print(classification_report(y_train,predicted_train))\n",
    "print()\n",
    "print('test set performance')\n",
    "print(classification_report(y_test,predicted_test))\n",
    "print()\n",
    "print('test set confusion matrix')\n",
    "cm = confusion_matrix(y_test,predicted_test)\n",
    "disp = ConfusionMatrixDisplay(cm)\n",
    "disp.plot();"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Well what if instead of train/test, we do cross validation?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.model_selection import cross_val_predict\n",
    "\n",
    "print('lets do 5 fold cross validation...')\n",
    "X_xval, X_test, y_xval, y_test = train_test_split( X, y, test_size=0.3, random_state=777)\n",
    "\n",
    "print('xval/test set sizes')\n",
    "print(X_xval.shape, X_test.shape, y_xval.shape, y_test.shape)\n",
    "print()\n",
    "\n",
    "pclassifier = KNeighborsClassifier(n_neighbors=15)\n",
    "\n",
    "predicted_xval = cross_val_predict(pclassifier, X_xval, y_xval, cv=5)\n",
    "predicted_test = pclassifier.fit(X_xval,y_xval).predict(X_test)\n",
    "\n",
    "print('xval performance')\n",
    "print(classification_report(y_xval,predicted_xval))\n",
    "print()\n",
    "print('test set performance')\n",
    "print(classification_report(y_test,predicted_test))\n",
    "print()\n",
    "print('test set confusion matrix')\n",
    "cm = confusion_matrix(y_test,predicted_test)\n",
    "disp = ConfusionMatrixDisplay(cm)\n",
    "disp.plot();"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Wow... but cross-validation is really useful for looking at the BEST set of hyperparameters possible... is there a way to do that?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.model_selection import GridSearchCV\n",
    "\n",
    "print('lets do grid search using 5 fold cross validation... and this time lets use all the measurements not just the best two')\n",
    "X = df[['bill_length_mm','bill_depth_mm','flipper_length_mm', 'body_mass_g']]\n",
    "y = df['species']\n",
    "\n",
    "preproc = StandardScaler()\n",
    "\n",
    "preproc.fit(X)\n",
    "scaled_data = preproc.transform(X)\n",
    "X['bill_length_mm'] = scaled_data[:,0]\n",
    "X['bill_depth_mm'] = scaled_data[:,1]\n",
    "X['flipper_length_mm'] = scaled_data[:,2]\n",
    "X['body_mass_g'] = scaled_data[:,3]\n",
    "\n",
    "X_xval, X_test, y_xval, y_test = train_test_split( X, y, test_size=0.3, random_state=777)\n",
    "\n",
    "print('xval/test set sizes')\n",
    "print(X_xval.shape, X_test.shape, y_xval.shape, y_test.shape)\n",
    "print()\n",
    "\n",
    "# note we don't need to put n_neightbors here now :)\n",
    "pclassifier = KNeighborsClassifier()\n",
    "params = { 'n_neighbors':[1,3,5,7,9,11,13,15]}\n",
    "\n",
    "searcher = GridSearchCV(pclassifier, param_grid=params, cv=5)\n",
    "\n",
    "# at this point we run a bunch of cross validations\n",
    "searcher.fit(X_xval, y_xval)\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# a fitted set of cvs, one for each param set in the grid\n",
    "searcher.cv_results_"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# you can match this up with mean_test_score above\n",
    "searcher.best_params_"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# searcher, when used as a classifier, uses the best_params_ by default\n",
    "predicted_test = searcher.fit(X_xval,y_xval).predict(X_test)\n",
    "\n",
    "\n",
    "print('test set performance')\n",
    "print(classification_report(y_test,predicted_test))\n",
    "print()\n",
    "print('test set confusion matrix')\n",
    "cm = confusion_matrix(y_test,predicted_test)\n",
    "disp = ConfusionMatrixDisplay(cm)\n",
    "disp.plot();"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Support Vector Machines: Maximizing the *Margin*\n",
    "\n",
    "Support vector machines are another ML technique\n",
    "\n",
    "The intuition is this: rather than simply drawing a zero-width line between the classes, we can draw around each line a *margin* of some width, up to the nearest point.\n",
    "Here is an example of how this might look:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.datasets import make_blobs\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "X, y = make_blobs(n_samples=50, centers=2,\n",
    "                  random_state=0, cluster_std=0.60)\n",
    "\n",
    "xfit = np.linspace(-1, 3.5)\n",
    "sns.scatterplot(x=X[:, 0], y=X[:, 1], hue=y, s=75);\n",
    "\n",
    "for m, b, d in [(1, 0.65, 0.33), (0.5, 1.6, 0.55), (-0.2, 2.9, 0.2)]:\n",
    "    yfit = m * xfit + b\n",
    "    plt.plot(xfit, yfit, '-k')\n",
    "    plt.fill_between(xfit, yfit - d, yfit + d, edgecolor='none',\n",
    "                     color='#AAAAAA', alpha=0.4)\n",
    "\n",
    "plt.xlim(-1, 3.5);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In support vector machines, the line that maximizes this margin is the one we will choose as the optimal model.\n",
    "Support vector machines are an example of such a *maximum margin* estimator."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Fitting a support vector machine\n",
    "\n",
    "Let's see the result of an actual fit to this data: we will use Scikit-Learn's support vector classifier to train an SVM model on this data.\n",
    "For the time being, we will use a linear kernel and set the ``C`` parameter to a very large number (we'll discuss the meaning of these in more depth momentarily)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.svm import SVC # \"Support vector classifier\"\n",
    "model = SVC(kernel='linear', C=1E10)\n",
    "model.fit(X, y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dir(model)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To better visualize what's happening here, let's create a quick convenience function that will plot SVM decision boundaries for us:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_svc_decision_function(model, ax=None, plot_support=True):\n",
    "    \"\"\"Plot the decision function for a 2D SVC\"\"\"\n",
    "    if ax is None:\n",
    "        ax = plt.gca()\n",
    "    xlim = ax.get_xlim()\n",
    "    ylim = ax.get_ylim()\n",
    "    \n",
    "    # create grid to evaluate model\n",
    "    x = np.linspace(xlim[0], xlim[1], 30)\n",
    "    y = np.linspace(ylim[0], ylim[1], 30)\n",
    "    Y, X = np.meshgrid(y, x)\n",
    "    xy = np.vstack([X.ravel(), Y.ravel()]).T\n",
    "    P = model.decision_function(xy).reshape(X.shape)\n",
    "    \n",
    "    # plot decision boundary and margins\n",
    "    ax.contour(X, Y, P, colors='k',\n",
    "               levels=[-1, 0, 1], alpha=0.5,\n",
    "               linestyles=['--', '-', '--'])\n",
    "    \n",
    "    # plot support vectors\n",
    "    if plot_support:\n",
    "        ax.scatter(model.support_vectors_[:, 0],\n",
    "                   model.support_vectors_[:, 1],\n",
    "                   s=300, linewidth=3, facecolors='none', edgecolor='k');\n",
    "    ax.set_xlim(xlim)\n",
    "    ax.set_ylim(ylim)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sns.scatterplot(x=X[:, 0], y=X[:, 1], hue=y, s=75);\n",
    "\n",
    "plot_svc_decision_function(model);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This is the dividing line that maximizes the margin between the two sets of points.\n",
    "Notice that a few of the training points just touch the margin: they are indicated by the black circles in this figure.\n",
    "These points are the pivotal elements of this fit, and are known as the *support vectors*, and give the algorithm its name.\n",
    "In Scikit-Learn, the identity of these points are stored in the ``support_vectors_`` attribute of the classifier:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model.support_vectors_"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A key to this classifier's success is that for the fit, only the position of the support vectors matter; any points further from the margin which are on the correct side do not modify the fit!\n",
    "Technically, this is because these points do not contribute to the loss function used to fit the model, so their position and number do not matter so long as they do not cross the margin.\n",
    "\n",
    "We can see this, for example, if we plot the model learned from the first 60 points and first 120 points of this dataset:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "\n",
    "def plot_svm(N=10, ax=None):\n",
    "    X, y = make_blobs(n_samples=200, centers=2,\n",
    "                      random_state=0, cluster_std=0.60)\n",
    "    X = X[:N]\n",
    "    y = y[:N]\n",
    "    model = SVC(kernel='linear', C=1E10)\n",
    "    model.fit(X, y)\n",
    "    \n",
    "    ax = ax or plt.gca()\n",
    "    sns.scatterplot(x=X[:, 0], y=X[:, 1], hue=y, s=75, ax=ax)\n",
    "    ax.set_xlim(-1, 4)\n",
    "    ax.set_ylim(-1, 6)\n",
    "    plot_svc_decision_function(model, ax)\n",
    "\n",
    "fig, ax = plt.subplots(1, 3, figsize=(16, 6))\n",
    "fig.subplots_adjust(left=0.0625, right=0.95, wspace=0.1)\n",
    "for axi, N in zip(ax, [30, 60, 120]):\n",
    "    plot_svm(N, axi)\n",
    "    axi.set_title('N = {0}'.format(N))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the left panel, we see the model and the support vectors for just 10 training points, in the middle we see  60 training points and on the right panel 120 points.  While going from 30 to 60 caused changes, going from 60 to 120 did not becuase the support vectors (the points on the margin line) had not changed.\n",
    "\n",
    "This insensitivity to the exact behavior of distant points is one of the strengths of the SVM model."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Beyond linear boundaries: Kernel SVM\n",
    "\n",
    "Where SVM becomes extremely powerful is when it is combined with *kernels*.\n",
    "We have seen a version of kernels before: Remember polynomial features from the regression lectures?  \n",
    "\n",
    "There we projected our data into higher-dimensional space defined by polynomials functions, and thereby were able to fit for nonlinear relationships with a linear regression.\n",
    "\n",
    "In SVM models, we can use a version of the same idea.\n",
    "To motivate the need for kernels, let's look at some data that is not linearly separable:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.datasets import make_circles\n",
    "X, y = make_circles(100, factor=.1, noise=.1)\n",
    "\n",
    "clf = SVC(kernel='linear').fit(X, y)\n",
    "\n",
    "sns.scatterplot(x=X[:, 0], y=X[:, 1], hue=y, s=75);\n",
    "plt.xlabel('x1')\n",
    "plt.ylabel('x2');\n",
    "#plot_svc_decision_function(clf, plot_support=False);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It is clear that no linear discrimination will *ever* be able to separate this data.\n",
    "\n",
    "One simple projection we could use would be to compute a *radial basis function* centered on the middle clump:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "phi = np.exp(-(X ** 2).sum(1))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Adding phi as a 3rd dimension gives us this:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib notebook\n",
    "\n",
    "from mpl_toolkits import mplot3d\n",
    "from matplotlib.colors import ListedColormap\n",
    "cmap = ListedColormap(sns.color_palette()[:2])\n",
    "\n",
    "ax = plt.subplot(projection='3d')\n",
    "ax.scatter3D(X[:, 0], X[:, 1], phi, c=y, s=7, cmap=cmap)\n",
    "ax.set_xlabel('x1')\n",
    "ax.set_ylabel('x2')\n",
    "ax.set_zlabel('$\\phi(\\mathbf{x})_3$');\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can see that with this additional dimension, the data becomes trivially linearly separable, by drawing a separating plane at, say, *phi*=0.7.\n",
    "\n",
    "Here we had to choose and carefully tune our projection: if we had not centered our radial basis function in the right location, we would not have seen such clean, linearly separable results.\n",
    "In general, the need to make such a choice is a problem: we would like to somehow automatically find the best basis functions to use.\n",
    "\n",
    "One strategy to this end is to compute a basis function centered at *every* point in the dataset, and let the SVM algorithm sift through the results.\n",
    "This type of basis function transformation is known as a *kernel transformation*, as it is based on a similarity relationship (or kernel) between each pair of points.\n",
    "\n",
    "A potential problem with this strategy—projecting $N$ points into $N$ dimensions—is that it might become very computationally intensive as $N$ grows large.\n",
    "However, because of a neat little procedure known as the [*kernel trick*](https://en.wikipedia.org/wiki/Kernel_trick), a fit on kernel-transformed data can be done implicitly—that is, without ever building the full $N$-dimensional representation of the kernel projection!\n",
    "This kernel trick is built into the SVM, and is one of the reasons the method is so powerful.\n",
    "\n",
    "In Scikit-Learn, we can apply kernelized SVM simply by changing our linear kernel to an RBF (radial basis function) kernel, using the ``kernel`` model hyperparameter:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "clf = SVC(kernel='poly', C=1)\n",
    "clf.fit(X, y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sns.scatterplot(x=X[:, 0], y=X[:, 1], hue=y, s=75)\n",
    "plot_svc_decision_function(clf)\n",
    "plt.scatter(clf.support_vectors_[:, 0], clf.support_vectors_[:, 1],\n",
    "            s=300, lw=1, facecolors='none', edgecolor='k');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using this kernelized support vector machine, we learn a suitable nonlinear decision boundary.\n",
    "This kernel transformation strategy is used often in machine learning to turn fast linear methods into fast nonlinear methods, especially for models in which the kernel trick can be used."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Tuning the SVM: Softening Margins\n",
    "\n",
    "Our discussion thus far has centered around very clean datasets, in which a perfect decision boundary exists.\n",
    "But what if your data has some amount of overlap?\n",
    "For example, you may have data like this:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X, y = make_blobs(n_samples=100, centers=2,\n",
    "                  random_state=0, cluster_std=1.2)\n",
    "sns.scatterplot(x=X[:, 0], y=X[:, 1], hue=y, s=75);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To handle this case, the SVM implementation has a bit of a fudge-factor which \"softens\" the margin: that is, it allows some of the points to creep into the margin if that allows a better fit.\n",
    "The hardness of the margin is controlled by a tuning parameter, most often known as $C$.\n",
    "For very large $C$, the margin is hard, and points cannot lie in it.\n",
    "For smaller $C$, the margin is softer, and can grow to encompass some points.\n",
    "\n",
    "The plot shown below gives a visual picture of how a changing $C$ parameter affects the final fit, via the softening of the margin:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X, y = make_blobs(n_samples=100, centers=2,\n",
    "                  random_state=0, cluster_std=0.8)\n",
    "\n",
    "fig, ax = plt.subplots(1, 2, figsize=(16, 6))\n",
    "fig.subplots_adjust(left=0.0625, right=0.95, wspace=0.1)\n",
    "\n",
    "for axi, C in zip(ax, [10.0, 0.1]):\n",
    "    model = SVC(kernel='linear', C=C).fit(X, y)\n",
    "    sns.scatterplot(x=X[:, 0], y=X[:, 1], hue=y, s=75, ax=axi)\n",
    "    plot_svc_decision_function(model, axi)\n",
    "    axi.scatter(model.support_vectors_[:, 0],\n",
    "                model.support_vectors_[:, 1],\n",
    "                s=300, lw=1, facecolors='none', edgecolor='k');\n",
    "    axi.set_title('C = {0:.1f}'.format(C), size=14)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The optimal value of the $C$ parameter will depend on your dataset, and should be tuned using cross-validation or a similar model selection procedure.\n",
    "\n",
    "\n",
    "In fact in general we would want to do model selection on all the hyperparameters of the SVM... $C$, what kernel to use, and whatever hyperparameters are unique to each kernel like polynomial order in $\\phi$.\n",
    "\n",
    "In fact sklearn gives us a powerful tool to do exactly this. Which I will demonstrate below.\n",
    "\n",
    "If you're actually following us all the way to this point you are getting into the territory of where I want you to be during this point in the quarter in COGS118A.  So if you are NOT following don't worry too much about it... this is advanced stuff "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# NOTE: this will take around 2-3 minutes to run\n",
    "\n",
    "from sklearn.compose import ColumnTransformer\n",
    "from sklearn.compose import make_column_selector as selector\n",
    "from sklearn.preprocessing import OneHotEncoder, StandardScaler\n",
    "from sklearn.pipeline import Pipeline\n",
    "from sklearn.svm import SVC\n",
    "from sklearn.model_selection import train_test_split, GridSearchCV\n",
    "\n",
    "\n",
    "# what we want to predict\n",
    "y = df['species']\n",
    "# the data we will predict from... note that we will include SEX\n",
    "# because penguin sex is actually pretty important in their size measurements\n",
    "X = df.drop(columns=['island','species'])\n",
    "\n",
    "# sex is a categorical variable... it needs to be one hot encoded!\n",
    "# we know from before that we need to z-score the numeric variables\n",
    "\n",
    "# here we setup a way to transform the categorical variables into one hot, \n",
    "# and to z-score the numeric variables.... completely automatically!!\n",
    "# remember how we had to do this with a bunch of different\n",
    "# steps manually the last time around??\n",
    "# if you'd like to understand how to do this more generally\n",
    "# search the sklearn docs for ColumnTransformer!\n",
    "numerical_columns_selector = selector(dtype_exclude=object)\n",
    "categorical_columns_selector = selector(dtype_include=object)\n",
    "\n",
    "numerical_columns = numerical_columns_selector(X)\n",
    "categorical_columns = categorical_columns_selector(X)\n",
    "\n",
    "onehot = OneHotEncoder()\n",
    "scaler = StandardScaler()\n",
    "\n",
    "preprocessor = ColumnTransformer([\n",
    "    ('onehot', onehot, categorical_columns),\n",
    "    ('zscore', scaler, numerical_columns)])\n",
    "\n",
    "\n",
    "# Create a pipeline \n",
    "# a pipeline puts together two (or 3 or however many) steps into one sklearn\n",
    "# estimator. That means we can use pipe instead of svc in any xval or grid search\n",
    "# and get the preprocessing done for us in the same step\n",
    "pipe = Pipeline([('make_features', preprocessor),\n",
    "                 ('classifier', SVC())])\n",
    "\n",
    "# SVM have very different hyperparameters depending on which kernel\n",
    "# potential hyperparams are gamma, coef0, degree and kernel.\n",
    "# here <x,x'> and ||x-x'|| are similarities between two vectors x,x'\n",
    "# linear kernel:     <x,x'>  [no hyperparams]\n",
    "# poly kernel:       ( gamma * <x,x'> + coef0 )^degree \n",
    "# rbf kernel:        exp( - gamma * ||x-x'||^2 )\n",
    "# sigmoid kernel:    tanh( gamma * <x,x'> + coef0 )\n",
    "\n",
    "# so how can we search using only using hypers appropriate for each kernel?\n",
    "# like this...make a list of different search spaces you want to execute\n",
    "search_space = [{'classifier__kernel': ['linear'],\n",
    "                 'classifier__C': np.logspace(-3, 2, 11) #11 steps between 10^-3 to 10^2\n",
    "                },\n",
    "                {'classifier__kernel': ['poly'],\n",
    "                 'classifier__gamma': np.logspace(-3, 2, 11),\n",
    "                 'classifier__degree': range(2,14),\n",
    "                 # NOPE... this isnt often important 'classifier__coef0': np.logspace(-3, 2, 11),\n",
    "                 'classifier__C': np.logspace(-3, 2, 11)\n",
    "                },\n",
    "                {'classifier__kernel': ['rbf'],\n",
    "                 'classifier__gamma': np.logspace(-3, 2, 11),\n",
    "                 'classifier__C': np.logspace(-3, 2, 11)\n",
    "                },\n",
    "                ]\n",
    "\n",
    "\n",
    "\n",
    "# Create a grid search object to find the best model\n",
    "best_model = GridSearchCV(pipe, search_space, cv=5, verbose=1)\n",
    "# play with different verbose=??\n",
    "\n",
    "# now lets split off some of our data for testing set... the rest will be used from cross validation\n",
    "X_train, X_test, y_train, y_test = train_test_split(X, y, \n",
    "                                                    test_size=0.25, random_state=101)\n",
    "\n",
    "\n",
    "# Fit grid search\n",
    "%time best_model.fit(X_train, y_train)\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# lets find out how well we do on the test set...\n",
    "print(best_model.best_params_)\n",
    "yhat = best_model.predict(X_test)\n",
    "print(classification_report(y_test, yhat))\n",
    "ConfusionMatrixDisplay.from_predictions(y_test, yhat);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "name": ""
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}