Bachelor Starter Code.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Hello, world!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Welcome to the Data Science Club's first prediction contest of the year! The club board wrote this starter code to give beginners to programming some already working code examples as a foundation for this contest.\n",
    "\n",
    "Our goals for this project:\n",
    " - Expose club members to important data science concepts and terminology. Since this field is so diverse, rapidly changing and growing, there exists a lot of overlapping terminology with many terms meaning the same thing (i.e. predictor variable vs. input variable vs. feature vs. independent variable). We've done our best to use consistent terminology to avoid confusion in this notebook!\n",
    " - Provide working code that serves as a foundation for future projects/exploration of this project\n",
    " - Take concepts from data science and apply it to something fun and outside the box! Who knew that predicting real love on reality TV would be so difficult?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src = \"images/machine_learning_xkcd.png\">"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "While to a lot of us data science feels like trying to understand amortization schedules in FY accounting, we hope to help demystify some of these concepts using intuitive and relevant examples."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Table of Contents:\n",
    "* [The Bachelor Forecasting Challenge](#The-Bachelor-Forecasting-Challenge)\n",
    " * [Goals of Statistical Learning](#Goal-of-Statistical-Learning)\n",
    " * [Prediction and Inference](#Prediction-and-Inference)\n",
    "* [Python](#Python)\n",
    " * [Load Packages](#Load-Packages)\n",
    " * [Load Data](#Load-Data)\n",
    "* [Data Cleaning and Manipulation](#Data-Cleaning-and-Manipulation)\n",
    "* [Exploratory Data Analysis](#Exploratory-Data-Analysis)\n",
    "* [Linear Regression](#Linear-Regression)\n",
    " * [R-squared and Adjusted R-squared](#R-squared-&-Adjusted-R-squared)\n",
    " * [Regression Coefficients](#Regression-Coefficients)\n",
    " * [Standard Error and Hypothesis Tests](#Standard-Error-and-Hypothesis-Tests)\n",
    " * [p-values](#p-values )\n",
    " * [Interpretation of Model Summary](#Interpretation-of-Model-Summary)\n",
    "* [Penalized Regressions](#Penalized-Regressions)\n",
    " * [Model Complexity and the Bias-Variance Tradeoff](#Model-Complexity-and-the-Bias-Variance-Tradeoff)\n",
    " * [Penalization / Regularization](#Penalization-/-Regularization)\n",
    " * [Standardized Data](#Standardized-Data)\n",
    " * [Model Parameters vs. Complexity Parameters](#Model-Parameters-vs.-Complexity-Parameters)\n",
    " * [Training, Validation, and Test Data](#Training,-Validation,-and-Test-Data)\n",
    " * [Shortcomings of Train/Validate Splits and Cross-Validation](#Shortcomings-of-Train/Validate-Splits-and-Cross-Validation)\n",
    " * [Lasso Regression](#Lasso-Regression)\n",
    " * [Summary of Penalized Regression](#Summary-of-Penalized-Regression)\n",
    "* [Ordinal Logistic Regression](#Ordinal-Logistic-Regression)\n",
    "* [Wrap Up](#Wrap-Up)\n",
    " * [Data Prep to Make Predictions](#Data-Prep-to-Make-Predictions)\n",
    " * [Making Predictions](#Making-Predictions)\n",
    " * [Save and Submit your Predictions](#Save-and-Submit-your-Predictions)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# The Bachelor Forecasting Challenge\n",
    "\n",
    "For our first project, we will be working with data found on contestants from ABC's [The Bachelor](https://en.wikipedia.org/wiki/The_Bachelor_(American_TV_series)). For this project, the data was obtained from [FiveThirtyEight's GitHub page](https://github.com/fivethirtyeight/data/tree/master/bachelorette) and from [data.world](https://data.world/aerispaha/the-bachelor-contestants). \n",
    "\n",
    "## Goal of Statistical Learning\n",
    "\n",
    "Our goal is to accurately predict the week that contestants are eliminated from the show (as early as Week 1 all the way to Week 10). We will be using contestant information such hometown, occupation, age, and performance on the first episode of each season to try and predict how long they will last on the show, otherwise known at their *Elimination Week*.\n",
    "\n",
    "In other words, our goal is to come up with a function that takes in observations of input/predictor variables (such as age, occupation, hometown) and outputs a prediction for a variable of interest (Elimination Week). \n",
    "\n",
    "Before diving in to the code and data, let's take a look at an example to motivate this problem!\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Example\n",
    " \n",
    "To use a B-school example, let's say that we want to forecast the sales of a product (our response variable) based upon our ad spend for that product and the ad spend our competitor's equivalent product (our two predictor variables). Here's some sample data for 5 months of sales:\n",
    "\n",
    "| Sales - $ Y $ | Our Ad Spend - $ X_1$ | Competitor's Ad Spend - $X_2$ |\n",
    "| :-: | :-: | :-: |\n",
    "| 62,543 | \\$1,045,000 | \\$1,932,000 |\n",
    "| 64,567 | \\$1,206,000 | \\$2,130,000 |\n",
    "| 59,905 | \\$1,113,000 | \\$2,110,000 |\n",
    "| 56,012 | \\$1,099,000 | \\$2,455,000 | \n",
    "| 63,996 | \\$1,305,000 | \\$2,005,000 |\n",
    "\n",
    "Using the above data, we want to estimate future sales as accurately as possible using the data we have on hand (our ad spend and our competitor's ad spend) we so can optimize our ad spend to maximize profits.\n",
    "\n",
    "Mathematically, this is represented as:\n",
    "\n",
    "$ \\hat{Y} = \\hat{f}(X) $\n",
    "\n",
    "where:\n",
    " - $ \\hat{Y} $ : prediction for our output variable (AKA **response variable**), so in this example, the response variable is unit sales\n",
    " - $ X $ : a matrix of our **predictor variables** (each variable is column in this table, each row is a month of data on ad spend)\n",
    " - $ \\hat{f}(X) $ : our **function** that will transform our input variables into our output variable\n",
    " \n",
    "**Note**: The little hats over the letters in the notation mean \"estimate\". We are estimating future values of $Y$ (our predictions) by coming up with our best estimate of the \"true relationship\" ($f$) between our predictor variables and our response variable. We denote this \"best estimate\" with $\\hat{f}$, since we are trying to estimate the true function itself. Think of it like when Isaac Newton was coming up with his best approximation for the law of gravitation. He found a nearly perfect representation of the gravitational attaction between two objects using their masses and distance between them as predictor variables. \n",
    "\n",
    "$ Y = \\hat{f}(m_1, m_2, R) $\n",
    "\n",
    "where\n",
    " - $Y$ is the gravitiational attraction between two objects\n",
    " - $m_1$ and $m_2$ are the masses of the two objects\n",
    " - $R$ is the distance between the two objects\n",
    " - $\\hat{f}$ is the estimate of the true relationship between these three predictor variables and the response variable\n",
    "  - $\\hat{f} = 9.81 \\frac{m_1 m_2}{R^2}$\n",
    "\n",
    "It wasn't until centuries later that Einstein found something that more was closer to the \"true relationship\" with his theory of general relativity.\n",
    "\n",
    "## Prediction and Inference\n",
    "\n",
    "A large part of machine learning is using different approaches to estimate the \"true function\" $f$ using past data. Why do we care about this? Estimating $f$ allows us to do things: **prediction** and **inference**.\n",
    "\n",
    "**Prediction** - We care less about interpretting the practical meaning of $\\hat{f}$ and care solely about the predictive accuracy of our model. \n",
    "\n",
    "**Inference** - We care more about interpretting the relationships between our predictor variables and our response variable. We could use this information to ask future strategic questions to improve outcomes, such as understanding the underlying relationship between ad spend and unit sales in order to maximize profits.\n",
    "\n",
    "Inference often goes hand in hand with **explainability**, or the relative ease of which a data scientist can explain how the model is making predictions and interpret the relationships between predictor variables and the response variable. Certain modeling methods are inherently more complex than other, but often can lead to more accurate results. Thus, we have this inverse relationship between accuracy and explainability, as shown below:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src = \"images/Accuracy_v_Explainability.jpg\">\n",
    "\n",
    "<h4 align=\"center\"> Source: Preet Gandhi - KDnuggets - https://www.kdnuggets.com/2019/01/explainable-ai.html </h4> "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For our notebook today, we will be looking at a few different **Regression Algorithms**:\n",
    " - Ordinary Least Squares (OLS) regression (AKA linear regression)\n",
    " - Penalized Regression (Lasso in particular)\n",
    " - Ordinal Regression\n",
    " \n",
    "These are all methods used to come up with an estimate $\\hat{f}$ for the true function $f$ using past data to \"train\" our model. Since we have a variety of different methods to try, we need a way to compare the performance of these methods against one another. To do so, we will be assessing the performance of our models using Mean Squared Error:\n",
    "\n",
    "$ MSE = \\frac{1}{n} \\sum_{i=1}^n (y_i - \\hat{y}_i)^2 $\n",
    "\n",
    "where:\n",
    " - $ n $ : the number of **observations**, or number of rows in our data\n",
    " - $ y $: the true value of **response variable**\n",
    " - $ \\hat{y} $: our **predicted value** of the response variable\n",
    " \n",
    "Without further ado, let's dive into the coding and the data!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Python"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In recent years, Python has become the primary language of choice for data scientists. This is largely due to:\n",
    " - The relative ease of learning the language compared to other programming languages\n",
    " - It's [open source](https://opensource.com/resources/what-open-source) so anyone can contribute to it and access it for free\n",
    " - Vast libraries of optimized code to make a data scientist's life easier\n",
    " \n",
    "<img src = 'images/python.png'>\n",
    "\n",
    "Another popular choice for data scientists is [R](https://www.r-project.org/about.html), but honestly I am using Python for this exercise because I don't wish that evil on anyone. Python is just way cooler. If you want a real analysis of the pros/cons of each langauge, check out this post on [KDnuggets](https://www.kdnuggets.com/2015/05/r-vs-python-data-science.html)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Load Packages\n",
    "\n",
    "Data scientists will usually start each notebook by importing a set of packages that contain already written code that will make our lives way easier. You can think of these packages like a collection of helper functions that an expert programmer already wrote so you don't have to! Each package that we import will serve a particular purpose. Some packages specialize in manipulating data (i.e. [pandas](https://pandas.pydata.org/) and [numpy](https://numpy.org/)), others specialize in visualizing data (i.e. [seaborn](https://seaborn.pydata.org/index.html) and [matplotlib](https://matplotlib.org/)), others specialize in machine learning applications ([sklearn](https://scikit-learn.org/stable/) and [tensorflow](https://www.tensorflow.org/)), etc. \n",
    "\n",
    "The one the main benefits of importing this code is so that we don't have to reinvent the wheel! Why write something like a function from scratch that someone else has already done (and probably done much better than we could do ourselves!). These packages are open source as well, so anyone can contribute to these to improve upon existing code or to add new functionality. \n",
    "\n",
    "These are the packages we will be using in this analysis. Data scientists typically import the packages of their choosing at the beginning of each notebook for organizational purposes. If they realize they need to import something new, they will scroll back up to the top and add a new import statement into this cell."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Note**: Not all packages are installed on this server, so you would be getting an error that says something like `Module not found`. You will first need to download and install the missing packages. Luckily, Python has a very easy way to do this! You can install any missing packages by running `!pip install <package_name>` in your notebook. Pip is a python package installer that will automatically find the source code for the package you want to use, download it, and install it."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Run this cell to install any missing packages, just replace \"pandas_profiling\" with the missing package name\n",
    "!pip install pandas_profiling geopy mord\n",
    "!pip install --upgrade seaborn"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#To manipulate datasets efficiently\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "\n",
    "#Do regular expression matching\n",
    "import re\n",
    "\n",
    "#Data visualization packages\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "\n",
    "#To do exploratory data analysis\n",
    "from pandas_profiling import ProfileReport\n",
    "\n",
    "#To do geospatial analysis\n",
    "import geopy\n",
    "from geopy.geocoders import Nominatim\n",
    "from geopy.extra.rate_limiter import RateLimiter\n",
    "from geopy import distance as gp_distance\n",
    "from geopy import Point\n",
    "\n",
    "#Regressions\n",
    "import statsmodels.api as sm\n",
    "from sklearn.linear_model import LinearRegression\n",
    "from sklearn.linear_model import Lasso\n",
    "from sklearn.linear_model import LassoCV\n",
    "import mord #for Ordinal Logistic Regression\n",
    "\n",
    "#More data manipulation and scoring metrics\n",
    "from sklearn.model_selection import train_test_split\n",
    "from sklearn.metrics import mean_absolute_error\n",
    "from sklearn.metrics import mean_squared_error\n",
    "from sklearn.metrics import zero_one_loss\n",
    "from sklearn.preprocessing import OneHotEncoder\n",
    "from sklearn.preprocessing import StandardScaler"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You'll notice two main types of syntax above:\n",
    " - `import <package> as <alias>`\n",
    " - `from <package> import <function>`\n",
    "\n",
    "The syntax `import <package> as <alias>` is telling python to load *all of the functions* that exist within a particular package. We can give the package an **alias**, or a shorthand name so we don't have to write `pandas` everytime we want to use a function from that package, we can simply write `pd`.\n",
    "\n",
    "The syntax `from <package> import <function>` is telling python to load *only one specific function* from the package. Some packages can contain hundreds of functions and we may not need all of them for our analysis. This allows us to have more of an a la carte approach than an all-you-can-eat approach."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Load Data\n",
    "\n",
    "### Contestants Data\n",
    "\n",
    "Now we'll start by loading in the data on our contestants from a csv (\"comma separated value\") file. Csv files are a common data format that data scientists work with, and thankfully `pandas` already has a built-in function that will read the csv file into a Pandas DataFrame.\n",
    "\n",
    "Pandas DataFrames are a tabular data structure that allow us to easily manipulate data in python. With pandas we can:\n",
    " - Load data in from various file formats (csv, Excel, txt)\n",
    " - Display our data\n",
    " - Create new columns in our dataset\n",
    " - Drop rows\n",
    " - Fill in missing values\n",
    " - Aggregate data into pivot tables\n",
    " - Join datasets together\n",
    " - Much, much more!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#This will load our csv file into a Pandas DataFrame\n",
    "contestants_df = pd.read_csv('bachelor-contestants.csv')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "When using a function like `read_csv` from the pandas package, we call the function with `<package_name/abbreviation>.<function_name>` or `pd.read_csv`. If we just imported one function from a package instead of the whole package, you can simply call the fuction with `read_csv`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Python functions often have a lot of different options built in to the function that accomplish certain tasks. You can easily read the documentation on a function by using a \"?\" before the function name. Some functions have more documentation than others, depending on how lazy the developers were :)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Get the documentation for a function\n",
    "?pd.read_csv"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "That was a lot of documentation! You won't always have to read this, but when working with new functions reading the documentation will help you understand how to use a function as the developers intended."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's look at the result here by printing out the first 5 rows using the `head` attribute from pandas."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "contestants_df.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can look at a specific number of rows by passing in the number of rows we wish to see into `head`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Look at the first 15 rows\n",
    "contestants_df.head(15)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The `summarize_dataframe` function will allow us to get a nice summary table of our data. This is an example of a user-defined function, originally written by Darden's Prof. Albert. I won't go into the nitty gritty in explaining the code here, but generally we define a function by giving it a name (`summarize_dataframe`) as well as the inputs it expects to recieve (just `df` in this case)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Define the function \n",
    "#Function originally written by Prof. Michael Albert\n",
    "def summarize_dataframe(df):\n",
    "    \"\"\"Summarize a dataframe, and report missing values.\"\"\"\n",
    "    missing_values = pd.DataFrame({'Variable Name': df.columns,\n",
    "                                   'Data Type': df.dtypes,\n",
    "                                   'Missing Values': df.isnull().sum(),\n",
    "                                   'Unique Values': [df[name].nunique() for name in df.columns]}\n",
    "                                 ).set_index('Variable Name')\n",
    "    with pd.option_context(\"display.max_rows\", 1000):\n",
    "        display(pd.concat([missing_values, df.describe(include='all').transpose()], axis=1).fillna(\"\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we can call this function and pass in any dataframe that we want to summarize! Let's use our `contestants_df`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "summarize_dataframe(contestants_df)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Some things to note:\n",
    " - The average age of a contestant on *The Bachelor* is 26 y/o\n",
    " - Average height is 5'5\", but the majority of our dataset is missing information on height\n",
    " - 50% of contestants are eliminated by week 2, and 75% eliminated by week 5"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The `ProfileReport` function from the `pandas_profiling` package gives us a great overview of our dataset without having to write too much code. Here we can see common values, the distributions of our data.\n",
    " "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "profile = ProfileReport(df = contestants_df, minimal = True\n",
    "                       )\n",
    "profile.to_notebook_iframe()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Other takeaways:\n",
    " - Age looks to be centered around 25, with contestants as young as 21 or as old as 36\n",
    " - The top professions of contestants are Attorneys (look out JD/MBAs!), Nannies, Flight Attendants, Teachers, and Students\n",
    " - Dallas, Chicago, Houston, LA, and Seattle are the top hometowns where contestants come from.\n",
    " - The histogram on seasons is misleading since some bins contain multiple seasons"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Bachelor Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "bachelor_df = pd.read_csv('bachelors.csv')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's get an overview of this data as well with our `summarize_dataframe` function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "summarize_dataframe(bachelor_df)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's also get a sense for the dimensions of our data, or the number of rows and columns. We can use the `.shape` attribute of a dataframe to print this information."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "bachelor_df.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can see that our bachelor dataframe has 16 rows and 5 columns."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since our dataset is so small, we can simply put `bachelor_df` to print the data out to the screen"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "bachelor_df"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Data Cleaning and Manipulation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In both datasets, height doesn't seem to be very consistently populated. For the sake of this analysis, let's drop this drop from each dataset"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Axis = 1 tells Pandas to look at column labels, Axis = 0 tells Pandas to look at row labels\n",
    "contestants_df = contestants_df.drop('Height', axis = 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "bachelor_df = bachelor_df.drop('Height', axis = 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "contestants_df"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's also try to calculate the distance between a contestant's and Bachelor's Hometowns to see if that has any affect on a constestant's success. This will take a few minutes to run, since it is calling an external [API](https://en.wikipedia.org/wiki/API) (just a way for two programs to talk to one another) for each row in the dataset.\n",
    "\n",
    "Don't worry about understanding the code here. Most of it was copy-pasted from [stack overflow](https://stackoverflow.com/questions) anyways. \n",
    "\n",
    "*Sidebar*: If you can't find an answer from looking at documentation, a good next step it to type your problem into Google. 99% of the time someone else has asked the same question you have on stack overflow, a programming forum, where expert programmers will (usually) kindly respond to your programming dilemma.\n",
    "\n",
    "<img src = 'images/stack_overflow.png'>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "locator = Nominatim(user_agent=\"BachelorGeocoder\")\n",
    "# 1 - conveneint function to delay between geocoding calls\n",
    "geocode = RateLimiter(locator.geocode, min_delay_seconds=0.01)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Run this for contestants\n",
    "\n",
    "# create location column\n",
    "contestants_df['location'] = contestants_df['Hometown'].apply(geocode)\n",
    "# create longitude, laatitude and altitude from location column (returns tuple)\n",
    "contestants_df['point'] = contestants_df['location'].apply(lambda loc: tuple(loc[1]) if loc else None)\n",
    "# split point column into latitude, longitude columns\n",
    "contestants_df[['latitude', 'longitude']] = pd.DataFrame(contestants_df['point'].tolist(), index=contestants_df.index)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Run this for Bachelors\n",
    "\n",
    "# create location column\n",
    "bachelor_df['location'] = bachelor_df['Hometown'].apply(geocode)\n",
    "# create longitude, laatitude from location column (returns tuple)\n",
    "bachelor_df['point'] = bachelor_df['location'].apply(lambda loc: tuple(loc[1]) if loc else None)\n",
    "# split point column into latitude, longitude columns\n",
    "bachelor_df[['latitude', 'longitude']] = pd.DataFrame(bachelor_df['point'].tolist(), index=bachelor_df.index)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "contestants_df = contestants_df.drop(['latitude','longitude','location'], axis = 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "bachelor_df = bachelor_df.drop(['latitude','longitude','location'], axis = 1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's see if whether or not a candidate is from the U.S. has an impact. We'll start by getting the state/country value from the contestant hometown"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "contestants_df['State/Country'] = contestants_df['Hometown'].str.split(',').str[1]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let's add some logic to flag contestants who are international:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Regex search terms for non-US countries\n",
    "non_us_terms = ['Quebec','Italy','British Columbia','Alberta','Pernambuco','England','Brazil','Ontario']\n",
    "    \n",
    "#Logic to check State/Country for regex terms. Returns 1 (US) or 0 (non=US)    \n",
    "contestants_df['International_Flag'] = ((contestants_df['State/Country'].str.contains('|'.join(non_us_terms))) | (contestants_df['State/Country'].isna())).astype(int)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "contestants_df = contestants_df.drop(['State/Country'], axis = 1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let's do the same for the Bachelors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "bachelor_df['State/Country'] = bachelor_df['Hometown'].str.split(',').str[1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Logic to check State/Country for regex terms. Returns 1 (US) or 0 (non=US)    \n",
    "bachelor_df['International_Flag'] = (bachelor_df['State/Country'].str.contains('|'.join(non_us_terms)) | bachelor_df['State/Country'].isna()).astype(int)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "bachelor_df = bachelor_df.drop(['State/Country'], axis = 1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's join these datasets together with `pd.merge` on season. But before doing so, let's rename some of these columns so that we don't have duplicates in our new dataframe."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Mapping of old names (on left) and new names (on right)\n",
    "contestant_rename_map = {'Name':'Contestant_Name',\n",
    "                        'Age':'Contestant_Age',\n",
    "                        'Hometown':'Contestant_Hometown',\n",
    "                        'point':'Contestant_Hometown_Coordinates',\n",
    "                        'International_Flag':'Contestant_International_Flag'}\n",
    "\n",
    "#Rename the columns in the dataframe, axis = 1 is telling Pandas to look at column labels\n",
    "contestants_df = contestants_df.rename(contestant_rename_map, axis = 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Mapping of old names (on left) and new names (on right)\n",
    "bachelor_rename_map = {'Name':'Bachelor_Name',\n",
    "                        'Age':'Bachelor_Age',\n",
    "                        'Hometown':'Bachelor_Hometown',\n",
    "                        'point':'Bachelor_Hometown_Coordinates',\n",
    "                        'International_Flag':'Bachelor_International_Flag'}\n",
    "\n",
    "#Rename the columns in the dataframe, axis = 1 is telling Pandas to look at column labels\n",
    "bachelor_df = bachelor_df.rename(bachelor_rename_map, axis = 1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let's merge these datasets together now that we have unique column names."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "merged_df = contestants_df.merge(bachelor_df, on = 'Season')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Looks like there is one contestant, Cosetta Blanca, whose age is null. Let's replace (AKA \"impute\") that null value with the median age of the other contestants."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "merged_df[merged_df['Contestant_Age'].isna()]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "merged_df.loc[merged_df['Contestant_Age'].isna(), 'Contestant_Age'] = merged_df['Contestant_Age'].median()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we have Bachelor and contestant data in one dataframe! Let's try adding some additional columns that may give us more information about potential relationships between a contestant and the Bachelor that season."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Calculate the age difference between a Bachelor and the contestant\n",
    "merged_df['Age_Difference'] = merged_df['Bachelor_Age'] - merged_df['Contestant_Age']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Exploratory Data Analysis"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's look at the distribution of age difference between the contestants and their respective bachelors. We can use the `seaborn` package with the alias `sns` to plot this data. We'll use the `histplot` function to plot a histogram of the age differences"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "g = sns.histplot(merged_df['Age_Difference'], binwidth = 1)\n",
    "g.set_xticks(range(int(merged_df['Age_Difference'].min()),int(merged_df['Age_Difference'].max())))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Above is a histogram showing the distribution of age differences between the Bachelor and the contestant. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's calculate the distance between contestants and the Bachelor for that season"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Calculate the distances between a contestant's hometown and the bachelor's\n",
    "merged_df['Hometown_Distance'] = merged_df.apply(lambda row: gp_distance.distance(row['Contestant_Hometown_Coordinates'], row['Bachelor_Hometown_Coordinates']).miles, axis = 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "g = sns.histplot(merged_df['Hometown_Distance'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Unsurprisingly, there are two modes to this distribution. The one on the left is (usually) associated with contestants from North America and the one on the right represents international contestants, with the exception of a couple of seasons where the Bachelor was from Europe:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "bachelor_df[bachelor_df['Bachelor_International_Flag']==1]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, let's create a column `Winner?` for the contestants if they win the Bachelor for that season"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "winner_lookup_table = merged_df.groupby(['Season'])['ElimWeek'].max().reset_index()\n",
    "\n",
    "#Drop the row for Season 11 since there were no winners on this season\n",
    "\n",
    "winner_lookup_table = winner_lookup_table[winner_lookup_table['Season']!=11]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "lookup_df = merged_df.merge(winner_lookup_table, how = 'left', on = 'Season').fillna(0)\n",
    "lookup_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Let's rename the columns to avoid confusion\n",
    "lookup_df = lookup_df.rename({'ElimWeek_x':'ElimWeek',\n",
    "                  'ElimWeek_y':'MaxElimWeek'}, axis = 1)\n",
    "\n",
    "#We'll cast this column as an integer using the astype function\n",
    "lookup_df['Winner?'] = (lookup_df['ElimWeek'] == lookup_df['MaxElimWeek']).astype(int)\n",
    "\n",
    "lookup_df = lookup_df.drop(['MaxElimWeek'], axis = 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "merged_df = lookup_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "merged_df"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let's look at aggregating our data. In pandas, we can use the `groupby` function to aggregate our data with a variety of aggregation functions (i.e. count, max, min, sum, average). \n",
    "\n",
    "Let's look at the categories of functions that contestants work in:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "merged_df.groupby(['Function'])['Function'].count()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let's get the number of times a contestant won a season by category:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "merged_df.groupby(['Function'])['Winner?'].sum()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here's the percetange of times contestants have won from each industry"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "merged_df.groupby(['Function'])['Winner?'].sum() / merged_df.groupby(['Function'])['Function'].count()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Right now, our `Function` column can't be interpretted by a regression algorithm since it is text and not numeric. We need to do something about that! How can we convert this information into a format that a regression algorithm would understand?\n",
    "\n",
    "We can convert this categorical data into numeric data by using something called \"dummy encoding\", or creating a column for each level in the category with a 1 if the contestant is in that category or 0 if they are not."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "len(merged_df['Function'].unique())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Our Function column has 15 distinct levels. So we will need to convert this single Function column into 15 separate columns with Yes/No indicators of if a contestant came from a particular industry. We will use the `OneHotEncoder` function from the `sklearn` package to easily to this for us, rather than having to write this logic ourselves."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Set up the dummy encoder\n",
    "dummy_encoder = OneHotEncoder()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Fit to our dataset (needs to be converted into a numpy array for this to work properly)\n",
    "dummy_encoder.fit(np.array(merged_df['Function']).reshape(-1,1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Creates our 15 distinct columns, one for each function, and does so for each contestant (423 contestants)\n",
    "sparse_matrix = dummy_encoder.transform(np.array(merged_df['Function']).reshape(-1,1)).toarray()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's double check the dimensions of this matrix to make sure we have the correct number of rows and columns:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sparse_matrix.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Sum across all rows to make sure we don't get any values != 1, since all contestants should only belong to 1 function in our dataset\n",
    "sparse_matrix.sum(axis = 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sparse_matrix"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Looks good! When working with new functions it's best practice to double check your assumptions with a few quick checks to make sure everything worked as intended."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Convert this into a pandas dataframe so we can join it back with the merged_df\n",
    "sparse_df = pd.DataFrame(data = sparse_matrix, columns = dummy_encoder.categories_)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "merged_df = merged_df.merge(sparse_df, left_index = True, right_index = True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "merged_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Clean up the column names by removing special characters like \"(\", \")\", \",\" , and \"''\"\n",
    "merged_df = merged_df.rename(columns = lambda x: re.sub(r'[(),\\']',r'',str(x)))\n",
    "\n",
    "#Replace spaces with underscores\n",
    "merged_df = merged_df.rename(columns = lambda x: re.sub(r' ',r'_',str(x)))\n",
    "\n",
    "#Replace slashes with underscores\n",
    "merged_df = merged_df.rename(columns = lambda x: re.sub(r'/',r'_',str(x)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's look at a correlation plot to visualize relationships between variables:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "corr = merged_df.corr()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mask = np.triu(np.ones_like(corr, dtype=bool))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f, ax = plt.subplots(figsize=(11, 9))\n",
    "cmap = sns.diverging_palette(230, 20, as_cmap=True)\n",
    "sns.heatmap(corr, mask=mask, cmap=cmap, vmax=1, vmin = -1, center=0,\n",
    "            square=True, linewidths=.5, cbar_kws={\"shrink\": .5})"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can also look at this just by printing the dataframe"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "corr.iloc[:10,:10]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can explore the relationships between pairs of variables by using a scatterplot. Seaborn has a useful function called `scatterplot` that will allow us to easily plot the relationship between two columns within a pandas dataframe."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sns.scatterplot(x = merged_df['Age_Difference'], y = merged_df['ElimWeek'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Seaborn also has a useful function called `regplot` that will not only plot the scatterplot of two variables, but will automatically calculate the line of best fit of the data and plot the confidence interval (parameter `ci`) around the line of best fit."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sns.regplot(x = merged_df['Hometown_Distance'], y = merged_df['ElimWeek'], ci = 95)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Linear Regression"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "After all that work, we finally have our data in a format where we can start building some models! Let's kick things off by looking at Ordinary Least Squares Regression, AKA Linear Regression.\n",
    "\n",
    "Linear Regressions fall under a broader class of models known as **Linear Models**. These models take the form of a weighted sum of predictor variables:\n",
    "\n",
    "$ y = \\beta_0 + \\sum \\beta_j x_j + \\epsilon $\n",
    "\n",
    "where:\n",
    " - $ x $ is a row of our $j$ predictor variables (i.e. a single contestants occupation, age, if they got a first impression rose or not, etc)\n",
    " - $ y $ is the response variable (the week a contestant is eliminated)\n",
    " - $ \\beta_0 $ is the intercept of the model, or the default value of the response variable when all predictor variables equal zero\n",
    " - $ \\beta_j $ is a coefficient for a predictor variable, with j total predictor variables, and indicates the relationship between a predictor variable and the response variable\n",
    " - $ \\epsilon $ is the error term of our model, which we assume is normally-distributed with a mean of zero\n",
    " \n",
    "To put this into terms of our model, let's say we want to use 3 predictor variables (AKA \"features\") to predict the elimination week:\n",
    " - Hometown Distance ($x_1$)\n",
    " - Age Difference ($x_2$)\n",
    " - First Impression Rose ($x_3$)\n",
    " \n",
    "Our linear model would look like this:\n",
    "\n",
    "$ y = \\beta_0 + \\beta_1 x_1 + \\beta_2 x_2 + \\beta_3 x_3 $\n",
    "\n",
    "The betas are the coefficients for each variable. Coefficients can be positive or negative and can us some information about the relationships between a predictor variable and the response variable.\n",
    "\n",
    "To find the optimal $ \\beta $ values, Linear Regression will pick $ \\beta $s that minimize the **residuals** of our model, or the distance between our model's predictions and the actual values. More specifically, Linear Regression will minimize the squared values of these residuals, or the **Residual Sum of Squares** (RSS):\n",
    "\n",
    "#### Math (totally optional)\n",
    "\n",
    "Residuals: $\\hat{\\epsilon}_i = y_i - \\hat{y}_i $\n",
    "\n",
    "where:\n",
    " - $i$ is a single observation in our training data, with $n$ total observations\n",
    " - $y$ is the actual value for our response variable in the training data\n",
    " - $\\hat{y}$ is our model's predction for the response variable given the training data\n",
    "\n",
    "Residual Sum of Squares: $ RSS = \\sum_{i}^n \\hat{\\epsilon}_i^2 $\n",
    "\n",
    "We can also rewrite the above equation to include our $\\beta$s, since our model's residuals will change depending on the $\\beta$ values chosen. This shows us that **RSS is minimized in linear regression by changing the $\\beta$s**:\n",
    "\n",
    "$ RSS(\\beta_0, \\beta_1, \\beta_2, \\beta_3) = \\sum_{i}^n (y_i - (\\beta_0 + \\beta_1 x_1 + \\beta_2 x_2 + \\beta_3 x_3) )^2 $\n",
    "\n",
    "#### Back to programming\n",
    "\n",
    "Let's get started using the `statmodels` package, which has some great functionality for creating linear regression models in python.\n",
    "\n",
    "We'll start by using the sub-package `OLS` from the statsmodels package, and from that we'll use the `from_forumla` function. This allows us to very easily define our linear regression model in terms of our response variable `ElimWeek` and our predictor variables `Hometown_Distance`, `Age_Difference`, and `First_Impression_Rose`. Note that the response variable is on the left-hand side of the formula and the predictor variables are on the right with a `+` sign between each one."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Here we are simply defining the structure of our model\n",
    "LinReg_Model1 = sm.OLS.from_formula(formula = 'ElimWeek ~ Hometown_Distance + Age_Difference + First_Impression_Rose', data = merged_df)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Now we are telling statsmodels to fit the model to our data\n",
    "LinReg_Model1_Results = LinReg_Model1.fit()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's write a quick function to calculate our model's Mean Squared Error:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def MSE(y_true, y_pred):\n",
    "    \"\"\"Returns the Mean Squared Error by passing the true values of the response variable (y_true) and a model's predicted values (y_pred)\"\"\"\n",
    "    if len(y_true) != len(y_pred):\n",
    "        print(\"Number of rows in predicted and true values do not match!\")\n",
    "    else:\n",
    "        n = len(y_pred)\n",
    "        MSE = ((y_true - y_pred)**2).sum() / n\n",
    "        print(MSE)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Generate our predicted values\n",
    "LinReg_Model1_Predictions = LinReg_Model1_Results.predict(merged_df)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Save our actual values into a variable\n",
    "Y_True = merged_df['ElimWeek']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "MSE(y_true = Y_True, y_pred = LinReg_Model1_Predictions)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Looks like our MSE is about 6.54. How does this compare with a very naive model, where we simply guess the average Elimation Week for every single contestant?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Naive_Predictions = np.repeat(merged_df['ElimWeek'].mean(), len(merged_df['ElimWeek']))\n",
    "\n",
    "MSE(y_true = merged_df['ElimWeek'], y_pred = Naive_Predictions)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So it looks like we're doing slightly better than guessing the average (3.4 weeks) for all contestants."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Now we can print and interpret the results with the summary function\n",
    "print(LinReg_Model1_Results.summary())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "That's a lot of information that the `summary` function gave us! Let's break this down step-by-step.\n",
    "\n",
    "## R-squared & Adjusted R-squared\n",
    "\n",
    "Let's start by looking at the `R-squared` and `Adj. R-squared`. In statistical terms, `R-squared` corresponds with how much variation in the data can be explained by the model we just built. In this case, about only 6.1% of the variability in the data can be explained by our model. The other 93.9% is explained by other information that we did not incorporate into our model, such as:\n",
    " - contestant height\n",
    " - if a contestant *is there for the right reasons*\n",
    " - about a billion other potential factors that might be able to explain \"love\" on a reality dating show\n",
    "\n",
    "`Adj. R-squared` AKA Adjusted R-squared is very similar to R-squared, but it penalizes the model for every new predictor variable that you add. We'll discuss why this is important later, but *generally we prefer to look at Adjusted R-squared rather than R-squared* so we can penalize more complex models and prevent overfitting our model to our data.\n",
    "\n",
    "## Regression Coefficients\n",
    "\n",
    "Next, let's turn our attention to the middle part of the summary where we see our Intercept term and 3 predictor variables. \n",
    "\n",
    "The `coef` column corresponds to the **estimates of the true coefficients** of each of our predictor variables in our model. In laymens terms, it means that a 1 unit change in a particular variable will correspond to a change in the response variable by the value of the cofficient. Using our model, obtaining the `First_Impression_Rose` would result in a contestant staying around 3.4495 weeks longer than if they had not gotten the rose in the first week.\n",
    "\n",
    "If we were to put this model into an equation, it would look like this:\n",
    "\n",
    "$ $ElimWeek$ = 3.4484 -0.0002($Hometown Distance$) + 0.0167($Age Difference$) + 3.4495($First Impression Rose$) $\n",
    "\n",
    "We could now use this formula to make predictions on our new contestants in the upcoming season!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's also calculate Mean Squared Error by defining our own function:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Standard Error and Hypothesis Tests\n",
    "\n",
    "The `std err` column tells us about the uncertainty in our estimate of each coefficient. Using the rule of thumb that $\\pm$ 2 standard deviations will cover about 95% of the possible estimates of each coefficient, we can infer a lot from looking at the standard error. Looking at the `First_Impression_Rose`, our model says that our estimate for the impact that getting the first impression rose has on a contestant's duration on the show could be $3.4495 \\pm 2*0.699$ weeks, or anywhere between 2.05 and 4.85 weeks. \n",
    "\n",
    "Looking at `Age_Difference`, our model says that the impact of a contestant having a 1 year larger age difference with the Bachlor could result in them staying between -0.0433 and 0.0767 weeks longer. This helps us build intuition that age difference probably does not make a big difference on how long a candidate stays on the show, since the model isn't even sure if they will stay longer or shorter because of an increase in age difference!\n",
    "\n",
    "What was just described above is a more intuitive way to think about **hypothesis tests**. In statistics, we use hypothesis tests to make assertions about relationships between variables. For example, we can write a **null hypothesis** (written as $H_0$) for each of our predictor variables in the model stating that *there is no relationship between the predictor and the response variables*. This is equivalent to saying that a regression coefficient $ \\beta_j = 0 $. We compare the null hypothesis with the **alternative hypothesis** (written as $H_A$), which would be the *there is some relationship between the predictor and response variable.*\n",
    "\n",
    "Mathematically, this is written as:\n",
    "\n",
    " - Null Hypothesis: $ H_0: \\beta_j = 0 $\n",
    " - Alt. Hypothesis: $ H_A: \\beta_j \\neq 0 $\n",
    " \n",
    "In other words, we are saying that some predictor variable has a relationship with the response variable or does not. So how do we quantify something like this? For that, we look to p-values."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## p-values\n",
    "\n",
    "The `P > |t|` is the **p-value** and tells us the probability that we would see a result as extreme or more extreme given a random sample. **Alternatively, it tells us the probability of us incorrectly rejecting our null hypothesis, which claims that there is no relationship between a given predictor variable and the response variable**. We can use the p-value to help us interpret which of our predictor variable(s) has a \"significant\" relationship with the response variable.\n",
    "\n",
    "Smaller p-values indicate that the chance of seeing such a relationship between the predictor and response variables due to random chance is incredibly unlikely, thus there must be a relationship of some kind between the two variables. When we see such a small p-value (common thresholds are below 0.05 or 0.01), we will reject the null hypothesis in favor of the alternative hypothesis - that there must be some kind of relationship between the two variables.\n",
    "\n",
    "**Note**: the choice of your *threshold for rejecting the null hypothesis* (also referred to as $ \\alpha $) is *subjective*. For some circumstances where is can be costly to make assumptions about relationships between a predictor and response variable, we would want to choose a very small threshold for significance. In other contexts where it is less costly, we might be OK using a higher theshold, such as 0.1. \n",
    "\n",
    "As a data scientist, you should pick your $ \\alpha $ before you run the model and interpret the results as to not bias your thinking."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Interpretation of Model Summary\n",
    "\n",
    "So, back to the interpretation of our results. We can see that getting the First Impression Rose has a strong impact on a contestant's ability to stay on the show from looking at the p-value for this predictor variable. We also see some evidence of the hometown distance contributing to a contestant's success on the show with a p-value of 0.093. While this is not less than 0.05, when dealing with challenging and complex problems like trying to predict love on a reality TV show, we might opt for using this as a feature to our model.\n",
    "\n",
    "Age difference, with a p-value of 0.584, does not seem to have a significant relationship Elimination Week, so let's try removing this variable from our model and see if it improves our Adjusted R-squared to something higher than 0.054."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Here we are simply defining the structure of our model, but remove Age_Difference as a predictor variable\n",
    "LinReg_Model2 = sm.OLS.from_formula(formula = 'ElimWeek ~ Hometown_Distance + First_Impression_Rose', data = merged_df)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Now we are telling statsmodels to fit the model to our data\n",
    "LinReg_Model2_Results = LinReg_Model2.fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "LinReg_Model2_Predictions = LinReg_Model2_Results.predict(merged_df)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Print out the MSE of the model\n",
    "MSE(Y_True, LinReg_Model2_Predictions)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Now we can print and interpret the results with the summary function\n",
    "print(LinReg_Model2_Results.summary())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Overall, the Adjusted R-squared of our model slighlty improved from 0.054 to 0.056. Looks like we still have a ways to go to improve the performance of our model!\n",
    "\n",
    "What if we decide to square the distance between a contestant's and Bachelor's hometowns? The idea behind making this transformation to this variable is that distance might be increasingly detrimental to a contestant the further the distance between their hometown and the Bachelor's due to cultural or regional differences.\n",
    "\n",
    "To build this into our model, we will use `numpy`'s `power` function to raise `Hometown_Distance` to the power of 2."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Here we are simply defining the structure of our model, but remove Age_Difference as a predictor variable\n",
    "LinReg_Model3 = sm.OLS.from_formula(formula = 'ElimWeek ~ np.power(Hometown_Distance, 2) + First_Impression_Rose', data = merged_df)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Now we are telling statsmodels to fit the model to our data\n",
    "LinReg_Model3_Results = LinReg_Model3.fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Now we can print and interpret the results with the summary function\n",
    "print(LinReg_Model3_Results.summary())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We see slight improvement by making this transformation to our predictor variable. The model's Adjusted R-squared increased to 0.058 and the p-value of our new feature is now 0.048. \n",
    "\n",
    "Finally, let's add the 15 dummy features that we built earlier that contain a contestant's job role. It wouldn't be an analysis by MBA students unless we had to ask where someone was working previously :)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Some logic to convert our categorical variables into a format that can be easily read into statmodel's formula format\n",
    "\n",
    "#Convert the categorical variables into a list\n",
    "categorical_var_list = dummy_encoder.categories_[0].tolist()\n",
    "\n",
    "#Replace the spaces with underscores\n",
    "categorical_var_list = [re.sub(r' ', r'_', str(var))  for var in categorical_var_list]\n",
    "\n",
    "#Replace slashes with underscores\n",
    "categorical_var_list = [re.sub(r'/', r'_', str(var))  for var in categorical_var_list]\n",
    "\n",
    "#Join elements of the list together into a single string with plus signs between each variable\n",
    "categorical_variables =  \" + \".join(categorical_var_list)\n",
    "\n",
    "categorical_variables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Here we are simply defining the structure of our model, now adding our categorical variables\n",
    "LinReg_Model4 = sm.OLS.from_formula(formula = 'ElimWeek ~ np.power(Hometown_Distance, 2) + First_Impression_Rose + ' + categorical_variables, data = merged_df)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Now we are telling statsmodels to fit the model to our data\n",
    "LinReg_Model4_Results = LinReg_Model4.fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Now we can print and interpret the results with the summary function\n",
    "print(LinReg_Model4_Results.summary())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It looks like many of these professions do not have significant relationships with the Elimination Week. Let's remove a bunch of these from the model and try one with only relatively small p-values (Corporate, Education, and Medical)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Here we are simply defining the structure of our model, now adding our categorical variables\n",
    "LinReg_Model5 = sm.OLS.from_formula(formula = 'ElimWeek ~ np.power(Hometown_Distance, 2) + First_Impression_Rose + Corporate + Education + Medical', data = merged_df)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Now we are telling statsmodels to fit the model to our data\n",
    "LinReg_Model5_Results = LinReg_Model5.fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Now we can print and interpret the results with the summary function\n",
    "print(LinReg_Model5_Results.summary())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So our Adjusted R-squared has increased slightly to 0.061 by removing the other features. Let's see if removing Medical, the feature with the highest p-value, will help use improve our Adjusted R-squared."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Here we are simply defining the structure of our model, now adding our categorical variables\n",
    "LinReg_Model6 = sm.OLS.from_formula(formula = 'ElimWeek ~ np.power(Hometown_Distance, 2) + First_Impression_Rose + Corporate + Education', data = merged_df)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Now we are telling statsmodels to fit the model to our data\n",
    "LinReg_Model6_Results = LinReg_Model6.fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "LinReg_Model6_Predictions = LinReg_Model6_Results.predict(merged_df)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Print out the MSE of the model\n",
    "MSE(Y_True, LinReg_Model6_Predictions)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Now we can print and interpret the results with the summary function\n",
    "print(LinReg_Model6_Results.summary())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Looks like removing the Medical variable lowered the Adjusted R-squared of the model. At this point, you could continue exploring which variables to include/remove to improve your model's Adjusted R-squared, but I will leave this to you to build upon further.\n",
    "\n",
    "## Summary of Linear Regression\n",
    "\n",
    "Pros:\n",
    " - Runs extremely fast on your data\n",
    " - Results are (relatively) easy to interpret\n",
    " \n",
    "Cons:\n",
    " - Inflexible compared to alternative methods, often underfits your data\n",
    " - Feature selection can be challenging, espeically when considering things like [p-hacking](https://en.wikipedia.org/wiki/Data_dredging)\n",
    " - Makes strong [assumptions](https://www.statisticssolutions.com/assumptions-of-linear-regression/) about the relationships between variables \n",
    " - Strongly impacted by outliers\n",
    " "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Penalized Regressions\n",
    "\n",
    "## Model Complexity and the Bias-Variance Tradeoff"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next up, we have another type of linear model, **Penalized Regression**. Penalized regression builds off of Linear Regression by including a *penalty term* for adding new features to our model. Why might we want to penalize ourselves for adding new features? Adding more features to our model increases the **complexity** of our model, which can be potentially detrimental to our modeling efforts as data scientists by overfitting our model to our data. Eliminating unnecessary complexity from our model helps us improve *prediction accuracy* on unseen data and *model inference* by eliminating unimportant features. Think of this concept like [Occam's Razor](https://en.wikipedia.org/wiki/Occam%27s_razor), where the simplest possible explanation is usually the best one.\n",
    "\n",
    "In other words, **our goal is to find the optimal level of complexity for our model to reduce model error**. Model error can be decomposed into three parts:\n",
    " - **Bias** - how far off a model's average prediction is from the actual value\n",
    " - **Varance** - how much a given prediction will vary depending on the training data we give our model\n",
    " - **Irredudible error** - error that we as data scientists can not feasibly control or take into account with our model (due to cost of collecting additional data, lack of understanding of the problem, etc)\n",
    " \n",
    "Mathematically, model error can be written as:\n",
    "\n",
    "$Model Error = Bias^2 + Variance + Irreducible Error$\n",
    "\n",
    "From looking at the above equation, we see that both bias and variance contribute to model error. As data scientists, we have the ability to control the bias and variance in our model and attempt to find the optimal tradeoff between the two.\n",
    "\n",
    "Here is a good way to visualize bias and variance:\n",
    "\n",
    "<img src = \"images/bias_variance_tradeoff.png\">\n",
    "<h4 align=\"center\"> Source: Scott Fortmann-Roe - http://scott.fortmann-roe.com/docs/BiasVariance.html </h4> \n",
    "\n",
    "\n",
    "If our model is too simple, it will underfit our data and have high **bias**. It will not be able to fully capture the necessary complexity of the problem due to it being overly simplistic. \n",
    "\n",
    "Conversely, if it's too complex, it will overfit to our data and have high **variance**. The model might see spurious relationships in the data that are actually meaningless in real life. It's kind like your one friend who looks at the label of a bottle of wine and all the sudden thinks he's a sommelier, then you give him a pour from a different bottle and he gives the same notes as before even though the first bottle was a red and this one is a white (we all have that friend).\n",
    "\n",
    "<img src = 'images/this_is_a_white.gif'>\n",
    "\n",
    "\n",
    "\n",
    "The tradeoff between bias and variance is creatively called the **bias-variance trade-off**. Underfit models tend to have high bias and low variance, whereas overfit models tend to how low bias and high variance. Ideally, we want to find the sweet spot by finding the optimal tradeoff between bias and variance.\n",
    "\n",
    "<img src = \"images/model_complexity.png\">\n",
    "<h4 align=\"center\"> Source: Scott Fortmann-Roe - http://scott.fortmann-roe.com/docs/BiasVariance.html </h4> \n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Penalization / Regularization\n",
    "\n",
    "Let's return our attention to penalized regression. So now we understand the importance of limiting model complexity in order to reduce variance in our model and prevent overfitting. Penalized regressions have a built-in way to limit complexity, namely by **shrinking** (AKA **regularizing**) the impact that a predictor variable has within a model by reducing the value of the coefficient. If you want a mathematic explanation, read below. If not, just remember that penalized regressions will shrink the values of our model's coefficients closer, if not all the way, to zero (depending on the type of shrinkage). (Insert bad George Costanza meme here).\n",
    "\n",
    "#### Math (totally optional)\n",
    "\n",
    "To contrast Linear Regression, what we just did above, with Penalized Regression, we can see how the two methods change how they select the coefficients or the $\\beta$s:\n",
    "\n",
    "*Basic Linear Regression* - We are minimizing Residual Sum of Squares (RSS) by selecting the optimal $\\beta$s:\n",
    "\n",
    "$ \\hat{\\beta} = \\underset{\\beta}{\\mathrm{argmin}} RSS(\\beta) $\n",
    "\n",
    "**Note**: \"$\\underset{\\beta}{\\mathrm{argmin}} $\" just means finding the $\\beta$s that minimize the RSS\n",
    "\n",
    "*Penalized Regression* - We are minimizing Residual Sum of Squares **plus a penalty term** by picking $\\beta$s and a **complexity parameter**  $\\lambda$ (AKA **tuning parameter**), which will control the strength of the penalty we apply for adding complexity to our model:\n",
    "\n",
    "$\\hat{\\beta} = \\underset{\\beta}{\\mathrm{argmin}} RSS(\\beta) + \\lambda P(\\beta) $\n",
    "\n",
    "where:\n",
    " - $\\lambda$ is the penalty tuning parameter, and will dictate how heavily we will penalize the model for adding new features/complexity\n",
    " - $P(\\beta)$ is a penalty function, of which there are number of kinds:\n",
    "   - Ridge Penalty: $ \\sum_{j=1}^p |{\\beta_j}|^2 $\n",
    "   - Lasso Penalty: $ \\sum_{j=1}^p {\\beta_j}^2 $\n",
    "\n",
    "### Summary\n",
    "\n",
    "Again, it's not so important to understand the math but to understand the more general concept of penalized regression. **We introduce the penalty term to automatically penalize our models for added complexity**. This will actually help our models better **generalize** to unseen data and find the optimal level of complexity for our model.\n",
    "\n",
    "Let's try implementing penalized regression in python!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Standardized Data\n",
    "\n",
    "To run our penalized regression techniques, we have to do a bit more data manipulation first. Because we are planning on shrinking some of the coefficients in our model to 0, we need to **standardize** our data so that all of our predictor variables are within roughly the same magnitude. For instance, `First_Impression_Rose` can only take on a value of 0 or 1, while `Hometown_Distance` can theoretically range from 0 to 1000s of miles. From our earlier example, we saw that `First_Impression_Rose` had a coefficient of about 3.5, while `Hometown_Distance` had a coefficient of -0.0002, yet both were significant at a level of $ \\alpha \\leq 0.1 $. The differences in magnitidue of the coefficients is largely because of the difference in scale between the two variables. Thus, it would be much easier to shrink `Hometown_Distance`'s coefficient to zero than `First_Impression_Rose`'s because it's so much closer to zero, and that's simply due to the difference in scale between the two predictor variables.\n",
    "\n",
    "To account for this, we need to standardize our data so each column has a mean of zero and a standard deviation of 1. That way, all of our data will have the same scale and we won't have to worry about the some variables being more easily penalized than others. To standardize our data, we subtract the mean and divide by the standard deviation for each predictor variable:\n",
    "\n",
    "standardized value $ = \\frac{x - \\mu}{\\sigma} $\n",
    "\n",
    "where:\n",
    " - $x$ is the data\n",
    " - $\\mu$ is the sample mean for a given predictor variable\n",
    " - $\\sigma$ is the sample standard deviation for a given predictor variable"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Here we define a list of continuous features that we want to scale\n",
    "#We do not need to standardize our binary features (i.e. First_Impression_Rose)\n",
    "Continuous_LinReg_features = ['Hometown_Distance','Age_Difference']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Define our scaling function from the sklearn package\n",
    "scaler = StandardScaler()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Fit our scaling function to the continuous features of our data\n",
    "scaler.fit(merged_df[Continuous_LinReg_features])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Create a dataset with our scaled data\n",
    "X_scaled = scaler.transform(merged_df[Continuous_LinReg_features])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Now we make a list of all of our categorical variables\n",
    "Binary_LinReg_features =  ['First_Impression_Rose'] + categorical_var_list"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Simply selecting the binary features of our data and converting into a numpy array\n",
    "X_binary = np.array(merged_df[Binary_LinReg_features])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Now we are merging our two datasets together into one\n",
    "X_penalized = np.concatenate((X_scaled, X_binary), axis = 1)\n",
    "Y = merged_df['ElimWeek']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Model Parameters vs. Complexity Parameters\n",
    "\n",
    "Now that we have our data scaled properly for penalized regression, we now want to split our data into **training** and **validation** sets. To understand what this means and why we do it, we need to understand the difference between **model parameters** and **complexity parameters**.\n",
    "\n",
    "Unlike with Linear Regression where we were only estimating **model parameters** (the $\\beta$ coefficients), we now also have to estimate the value of our **complexity parameter** $\\lambda$. We can think of complexity parameters a way to *change the amount of complexity in a model* whereas model parameters are ways that we optimally fit a model to our data *given a fixed level of complexity*. Complexity parameters are more high-level than model parameters and have a large impact on how well a model can generalize to new data. Because complexity parameters have such a large impact on model performance, we need to find the optimal value for our complexity parameter(s) *before* finding our optimal model parameters.\n",
    "\n",
    "## Training, Validation, and Test Data\n",
    "\n",
    "So how does any of this relate to splitting our data into training and validation sets? Remember, our goal is to build a model that performs best at predicting unseen data. Right now, we do not have any \"unseen data\" (AKA **test data**) since the new season has not finished and we don't know the outcomes of the show. Instead, we can *hold out* some of our data from our training set and keep it hidden from our models to simulate a test data set! We call the data we show our models **training data** and the data we withold from our models **validation data** since we are using this to validate or assess the performance of our different models and to find out which model performs best on unseen data.\n",
    "\n",
    " - **Training Data**: Data we show our model in order to find the best model parameters\n",
    " - **Validation Data**: Data we withold from our model to assess our models' performance and compare against one another\n",
    " - **Test Data**: Future data that we will feed into our best model to make predictions. This data will not have values for our response variable.\n",
    "\n",
    "Here's how this process generally works:\n",
    "\n",
    "1) Randomly split our data into training and validation sets\n",
    " - A rule of thumb is around 75-80% of your data should go into training and 20-25% should go into validation, but this can vary depending on your data\n",
    "\n",
    "2) Come up with a list of values that we want to test for $\\lambda$, our complexity parameter (`sklearn` will often do this for us!)\n",
    "\n",
    "3) Train our penalized regression models (one model per value of $\\lambda$) on our training data\n",
    " - The training set will be used to estimate the model parameters $\\beta$ in our penalized regression\n",
    " \n",
    "4) Compare the performance of each model against the validation set, which none of our models have seen\n",
    " - This will allow us to pick our optimal complexity parameter $\\lambda$. By assessing our models on data it has not seen, we can estimate how well each model **generalizes** to new data. This is exactly what we want to do when trying to predict the outcomes of new seasons.\n",
    "\n",
    "5) Finally refit the best model with the optimal complexity parameter to **both the training and validation data** so we can show the best model more data get even better estimates of our model parameters, $\\beta$.\n",
    "\n",
    "## Shortcomings of Train/Validate Splits and Cross-Validation\n",
    "\n",
    "Since we are randomly splitting our data, there's a chance that a random split could result in a biased split. For example, what if we randomly split our data and all of the winners of the Bachelor got put in the validation set? Our model would probably be much worse at predicting who would do well on the show since it would never get to train on any data where it saw a winner! How could we avoid this from happening? We can use something called **cross-validation**.\n",
    "\n",
    "Cross-validation addresses this problem by repeating the train/validate split multiple times to lower the chance of a \"bad split\" where we have imbalanced data between the training and validation sets. A common implementation of cross-validation is called K-fold cross validation, where K represents the number of times you repeat this process, and a \"fold\" represents $\\frac{1}{K}$ of your training data:\n",
    "\n",
    "<img src = 'images/kfold_cross_validation.png'>\n",
    "\n",
    "[Source](http://ethen8181.github.io/machine-learning/model_selection/model_selection.html) \n",
    "\n",
    "In the above example, K equals 5. So we first divide the training data into 5 equal groups. Next, we train our model in 4 of the 5 folds and then validate the performance of that model on the 5th fold. We repeat this process 4 more times, rotating which fold will be used as the validation set. Each model you test will be trained K times and thus will have K values of a performance metric (i.e. MSE). You can average each model's K MSE values to get an overall performance score for each model. \n",
    "\n",
    "The benefits of this method is that it ensures that each row of our data is used equally in training our model and that we do not end up with a biased validation set. The main drawback is that as you increase K, the time to assess your models' performance will greatly increase since you are training each model K times instead of once."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Split our data into training and test sets\n",
    "#Here we take 80% of our data into our training set and 20% goes into validation\n",
    "#Random state simply ensures that we get the same random split every time. This helps with reproducibility\n",
    "X_train, X_validate, Y_train, Y_validate = train_test_split(X_penalized, Y, test_size = 0.2, random_state = 37)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Lasso Regression\n",
    "\n",
    "Now that we have our data set up, our next step is to define our model. There are many types of penalized regression that we could use, but we are going to use a **Lasso Regression**, a specific type of Penalized Regression. Lasso has the benefit of completely reducing some of our variable's coefficients to zero. In other words, it tell us which variables we should keep and which ones we should get rid of!\n",
    "\n",
    "`sklean` also can do cross validation for us. The `LassoCV` function will create a Lasso Regression model AND implement cross validation at the same time! This makes our lives as data scientists much easier."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#cv: the number of folds we want to use in cross validation. 5 or 10 is usually a good choice.\n",
    "#Remember, this will help us find the optimal value of lambda\n",
    "LassoCV_model = LassoCV(cv = 10, normalize = False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Fit our model to our training data\n",
    "LassoCV_model.fit(X_train, Y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Prints out the coefficient from our regression (AKA the model parameters)\n",
    "LassoCV_coefficients = pd.DataFrame(zip(Continuous_LinReg_features + Binary_LinReg_features, LassoCV_model.coef_), columns = ['Feature','Coef'])\n",
    "LassoCV_coefficients"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Notice that many of the predictor variables now have coefficients of zero. Lasso has done some of the hard work for us in eliminating unimportant variables from our model!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "LassoCV_coefficients[LassoCV_coefficients['Coef']!=0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here's a visual representation of how our model's coefficent values change as we increase the value of our complexity parameter $\\lambda$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Code modified from Prof. Michael Albert's notebook on forecasting Clinton-Obama primary\n",
    "\n",
    "alphas = np.logspace(-5,0, 20)\n",
    "\n",
    "coefs = []\n",
    "for a in alphas:\n",
    "    lm_lasso = Lasso(alpha=a, fit_intercept=False, max_iter=1000000)\n",
    "    lm_lasso.fit(X_train, Y_train)\n",
    "    coefs.append(lm_lasso.coef_)\n",
    "    \n",
    "ax = plt.gca()\n",
    "\n",
    "for i in range(len(coefs[0])):\n",
    "    ax.plot(alphas, [coef[i] for coef in coefs], label=merged_df[Continuous_LinReg_features + Binary_LinReg_features].columns[i])\n",
    "ax.set_xscale('log')\n",
    "plt.xlabel('lambda')\n",
    "plt.ylabel('coefficient value')\n",
    "plt.title('Model parameters as a function of the strength of complexity parameter')\n",
    "plt.axis('tight')\n",
    "plt.legend(loc='right', bbox_to_anchor=(1.75, 0.5), shadow=True, ncol=1)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Notice that as our complexity parameter $\\lambda$ increases that many of our model parameters, the coefficents, start to shrink towards zero when $10^{-2} \\leq \\lambda \\leq 10^{-1}$. In other words, as we increase the penalty for complexity in our model, some of our predictor variables are effectively removed from our model, making it simpler.\n",
    "\n",
    "So what is our optimal value for $\\lambda$? Thankfully, `sklearn` will tell us this value automatically!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Prints out the optimal complexity parameter (sklearn calls it alpha but we've been referring to it as lambda)\n",
    "#It wouldn't be data science if no one could agree on terminology\n",
    "\n",
    "optimal_lambda = LassoCV_model.alpha_\n",
    "\n",
    "print(\"Optimal value for penalty term: \" +str(optimal_lambda))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now that we have our optimal value of our complexity parameter $\\lambda$, we have a good sense of the level of complexity that our model should have. Even better, `sklearn` will take the optimal value for our complexity parameter and refit our model to **all of our training data** to find the optimal estimates for our model parameters. So all we need to do now is generate predictions for our test set and estimate our model's performance on unseen data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Generate predictions for our model\n",
    "Lasso_model_predictions = LassoCV_model.predict(X_validate)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "MSE(Y_validate, Lasso_model_predictions)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Looks like Lasso does a slightly better job than linear regression! It's still not great, but honestly given how little data we have on the contestants, the Bachelors, and how challenging of a problem this is to model, it isn't a bad start!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Summary of Penalized Regression\n",
    "\n",
    "Penalized regression is a great next step from linear regression, as it allows us to better control the level of complexity in a model. It can even help us in automatically selecting which predictor variable we want to use in our model!\n",
    "\n",
    "Pros:\n",
    " - Can perform variable selection\n",
    " - Helps us control complexity, which prevents overfitting\n",
    " - Generally better predictive accuracy compared to Linear Regressions\n",
    "\n",
    "Cons:\n",
    " - Requires finding the complexity parameter\n",
    " - More effort to build a model / longer time to train model\n",
    " - Slightly harder to interpret / explain compared to a Linear Regression model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Ordinal Logistic Regression"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "One thing we've glossed over when reviewing our models is that our response variable, Elimination Week, is more of a ranking than a continuous number. However, the methods we've been using up until this point have been predicting non-integer values:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Lasso_model_predictions"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Obviously we can't say that a contestant will stay on a show for 3.12 weeks! They can only stay for 3 or 4 weeks. We could simply round these numbers to the nearest integer to compensate for this:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Lasso_model_predictions_rounded = np.round(Lasso_model_predictions)\n",
    "Lasso_model_predictions_rounded"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "MSE(Y_validate, Lasso_model_predictions_rounded)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "There is another way that we can do this with regression without rounding the result, we could use **Ordinal Logistic Regression**. In essence, Ordinal Logistic Regression is good for predicting ordinal response variables such as Movie Ratings, how an athlete might place in a swimming event, survey responses for on the [likert scale](https://en.wikipedia.org/wiki/Likert_scale), etc. In our case, we are equating a contestant's Elimination Week to a ranking in a competition, since in essence that's what a reality dating TV show is (sorry not sorry).\n",
    "\n",
    "Ordinal logistic regression looks at calculating the odds that an observation belongs to each ranking, and predicts the ranking that has the most likely value given the data. To use this kind of model, we need to make sure that our data passes the *proportional odds assumption*. To explain this, we first need to understand the difference between **probabilities** and **odds**. \n",
    "\n",
    "Probabilities take on values between 0 and 1 and represent the chance of an event occuring. \n",
    " - There will be a 40% chance of rain tomorrow or there is a 90% chance I get cold-called in DA next week.\n",
    "\n",
    "Odds can take on values between 0 and $\\infty$ and represent the ratio of an event occuring divided by an event not occuring. \n",
    " - E.g. I give it 3:1 odds that I messed something up on my DCF on my Finance final.\n",
    "\n",
    "#### Proportional Odds Assumption\n",
    "\n",
    "In order for us to use ordinal logistic regression, we need to make sure our data is well-suited for this kind of problem. Without getting too mathy, we want to make sure that our response variable, Elimination Week, is fairly evenly spread out. For the math, check out below (optional):\n",
    "\n",
    "#### Math (optional)\n",
    "\n",
    "We want to make sure the *logarithm of the odds* of making it a certain week should be equidistant, or evenly spaced from one another:\n",
    "\n",
    "Week 1: $log \\frac{p_1}{p_2 + p_3 + p_4 + ... + p_{11}} = X$\n",
    "\n",
    "Week 2: $log \\frac{p_1 + p_2}{p_3 + p_4 + p_5 + ... + p_{11}} = X + n$\n",
    "\n",
    "Week 3: $log \\frac{p_1 + p_2 + p_3}{p_4 + p_5 + p_6 ... + p_{11}} = X + 2n $\n",
    "\n",
    "...\n",
    "\n",
    "Week 10: $log \\frac{p_1 + p_2 + p_3 + p_{10}}{ p_{11}} X + 9n$\n",
    "\n",
    "Where $X$ can be any value and $n$ is a constant value that represents the equal distance between the log of the odds for each ranking. If $n$ is not consistent between rankings, then this is not a good choice of model for our data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Calculating the probabilities of getting eliminated by week by using our contestant data\n",
    "elim_prob = merged_df.groupby(['ElimWeek'])['ElimWeek'].count() / len(merged_df)\n",
    "elim_prob"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Create a list of values that represent the log odds for each ranking\n",
    "\n",
    "log_odds_list = []\n",
    "\n",
    "#Do the first ranking (week 1)\n",
    "log_odds_list.append(np.log(elim_prob.iloc[0].sum() / elim_prob.iloc[1:].sum()))\n",
    "\n",
    "#Do the other rankings (weeks 2-10)\n",
    "for i in range(1, 10, 1):\n",
    "    log_odds_list.append(np.log(elim_prob.iloc[0:i].sum() / elim_prob.iloc[i+1:].sum()))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Print the distances between log odds ratios to verify that n is fairly constant\n",
    "\n",
    "for i in range(0, len(log_odds_list)-1):\n",
    "    print(\"Difference between Weeks \" + str(i+1) + \" and \" +str(i+2) + \": \" + str(log_odds_list[i+1] - log_odds_list[i]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Excluding the difference between weeks 1 and 2, the other weeks are fairly consistent in terms of distance between the log odds ratios. Thus, we can proceed with modeling this problem as an Ordinal Logistic Regression.\n",
    "\n",
    "## Implementing Oridinal Logistic Regression\n",
    "\n",
    "We'll follow virtually the same steps as we did for the Linear Regression model, since the ordinal regression model we are using does not have any complexity parameters. We will be using the `mord` package for our Ordinal Logistic Regression model."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "features = ['Hometown_Distance','Age_Difference','First_Impression_Rose']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "features = features + categorical_var_list"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X = merged_df[features]\n",
    "Y = merged_df['ElimWeek']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_train, X_validate, Y_train, Y_validate = train_test_split(X, Y, test_size = 0.2, random_state = 37)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Ordinal_Logistic_Model = mord.LogisticAT()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Ordinal_Logistic_Model.fit(X_train, Y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "pd.DataFrame(zip(features, Ordinal_Logistic_Model.coef_), columns = ['Feature','Coef'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Ordinal_Logistic_Predictions = Ordinal_Logistic_Model.predict(X_validate)\n",
    "Ordinal_Logistic_Predictions"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that our model automatically produces integer outputs rather than fractions of a whole number."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "MSE(Y_validate, Ordinal_Logistic_Predictions)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Looks like this model does not perform as well as Penalized Regression given this test set. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Wrap Up"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Honestly that was way more material than we were planning to cover for this exercise! In this notebook, we covered:\n",
    "\n",
    " - The goals of statistical learning\n",
    " - How to use Python\n",
    " - Cleaning and manipulating data\n",
    " - Visualizing data\n",
    " - Building and interpretting Regression Models\n",
    " - Fundamental Data Science Concepts\n",
    "  - Bias-Variance Trade Off\n",
    "  - Model Complexity\n",
    "  - Assessing Model Performance\n",
    "  - Hypothesis Tests\n",
    "  - p-values\n",
    "  \n",
    "Give yourself a pat on the back for making it this far! Hopefully you're able to build an even better model that you can use to win your Bachelor league for years to come.\n",
    "\n",
    "For more resources and analytics on the Bachelor, be sure to check out some of these resources:\n",
    "\n",
    " - [@bachelordata](https://www.instagram.com/bachelordata/?hl=en)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Data Prep to Make Predictions\n",
    "\n",
    "So after all the hard work we've done, it's finally time to make predictions with our models! We first need to load in our **test data** for which we'll be making predictions on Season 25:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_df = pd.read_csv('test_data.csv')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "MattJames = pd.DataFrame({'Name':['Matt James'],\n",
    "                          'Age': [28],\n",
    "                          'Hometown': ['Raleigh, North Carolina'],\n",
    "                          'Season': [25]})\n",
    "\n",
    "MattJames"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We need to get our data into the same format as we used for our previous models. This means adding the same features (`Hometown_Distance`, `Age_Difference`, etc) and doing the same cleaning on this data as we did previously.\n",
    "\n",
    "**Note**: We just copy-pasted this code from above down here, replacing any mention of `contestant_df` with `test_df` and `bachelor_df` with `MattJames`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create location column\n",
    "test_df['location'] = test_df['Hometown'].apply(geocode)\n",
    "# create longitude, laatitude and altitude from location column (returns tuple)\n",
    "test_df['point'] = test_df['location'].apply(lambda loc: tuple(loc[1]) if loc else None)\n",
    "# split point column into latitude, longitude columns\n",
    "test_df[['latitude', 'longitude']] = pd.DataFrame(test_df['point'].tolist(), index=test_df.index)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Run this for Bachelors\n",
    "\n",
    "# create location column\n",
    "MattJames['location'] = MattJames['Hometown'].apply(geocode)\n",
    "# create longitude, laatitude from location column (returns tuple)\n",
    "MattJames['point'] = MattJames['location'].apply(lambda loc: tuple(loc[1]) if loc else None)\n",
    "# split point column into latitude, longitude columns\n",
    "MattJames[['latitude', 'longitude']] = pd.DataFrame(MattJames['point'].tolist(), index=MattJames.index)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_df = test_df.drop(['latitude','longitude','location'], axis = 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "MattJames = MattJames.drop(['latitude','longitude','location'], axis = 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_df['State/Country'] = test_df['Hometown'].str.split(',').str[1]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let's add some logic to flag contestants who are international:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Regex search terms for non-US countries\n",
    "non_us_terms = ['Quebec','Italy','British Columbia','Alberta','Pernambuco','England','Brazil','Ontario', 'Ethiopia']\n",
    "    \n",
    "#Logic to check State/Country for regex terms. Returns 1 (US) or 0 (non=US)    \n",
    "test_df['International_Flag'] = ((test_df['State/Country'].str.contains('|'.join(non_us_terms))) | (test_df['State/Country'].isna())).astype(int)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_df = test_df.drop(['State/Country'], axis = 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "MattJames['State/Country'] = MattJames['Hometown'].str.split(',').str[1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Logic to check State/Country for regex terms. Returns 1 (US) or 0 (non=US)    \n",
    "MattJames['International_Flag'] = (MattJames['State/Country'].str.contains('|'.join(non_us_terms)) | MattJames['State/Country'].isna()).astype(int)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "MattJames = MattJames.drop(['State/Country'], axis = 1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's join these datasets together with `pd.merge` on season. But before doing so, let's rename some of these columns so that we don't have duplicates in our new dataframe."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Mapping of old names (on left) and new names (on right)\n",
    "contestant_rename_map = {'Name':'Contestant_Name',\n",
    "                        'Age':'Contestant_Age',\n",
    "                        'Hometown':'Contestant_Hometown',\n",
    "                        'point':'Contestant_Hometown_Coordinates',\n",
    "                        'International_Flag':'Contestant_International_Flag'}\n",
    "\n",
    "#Rename the columns in the dataframe, axis = 1 is telling Pandas to look at column labels\n",
    "test_df = test_df.rename(contestant_rename_map, axis = 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Mapping of old names (on left) and new names (on right)\n",
    "bachelor_rename_map = {'Name':'Bachelor_Name',\n",
    "                        'Age':'Bachelor_Age',\n",
    "                        'Hometown':'Bachelor_Hometown',\n",
    "                        'point':'Bachelor_Hometown_Coordinates',\n",
    "                        'International_Flag':'Bachelor_International_Flag'}\n",
    "\n",
    "#Rename the columns in the dataframe, axis = 1 is telling Pandas to look at column labels\n",
    "MattJames = MattJames.rename(bachelor_rename_map, axis = 1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let's merge these datasets together now that we have unique column names."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "merged_test_df = test_df.merge(MattJames, on = 'Season')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Calculate the age difference between a Bachelor and the contestant\n",
    "merged_test_df['Age_Difference'] = merged_test_df['Bachelor_Age'] - merged_test_df['Contestant_Age']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "merged_test_df['Hometown_Distance'] = merged_test_df.apply(lambda row: gp_distance.distance(row['Contestant_Hometown_Coordinates'], row['Bachelor_Hometown_Coordinates']).miles, axis = 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Creates our 15 distinct columns, one for each function, and does so for each contestant (423 contestants)\n",
    "sparse_matrix = dummy_encoder.transform(np.array(merged_test_df['Function']).reshape(-1,1)).toarray()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Convert this into a pandas dataframe so we can join it back with the merged_df\n",
    "sparse_df = pd.DataFrame(data = sparse_matrix, columns = dummy_encoder.categories_)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "merged_test_df = merged_test_df.merge(sparse_df, left_index = True, right_index = True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Clean up the column names by removing special characters like \"(\", \")\", \",\" , and \"''\"\n",
    "merged_test_df = merged_test_df.rename(columns = lambda x: re.sub(r'[(),\\']',r'',str(x)))\n",
    "\n",
    "#Replace spaces with underscores\n",
    "merged_test_df = merged_test_df.rename(columns = lambda x: re.sub(r' ',r'_',str(x)))\n",
    "\n",
    "#Replace slashes with underscores\n",
    "merged_test_df = merged_test_df.rename(columns = lambda x: re.sub(r'/',r'_',str(x)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Making Predictions"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can make predictions with Linear Regression, Penalized Regression, and Ordinal Logistic Regression using the `predict` function."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Linear Regression predictions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "LinRegModel_Season25_Predictions = LinReg_Model6_Results.predict(merged_test_df)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Lasso Penalized Regression predictions (need to standardize our data first)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Here we define a list of continuous features that we want to scale\n",
    "#We do not need to standardize our binary features (i.e. First_Impression_Rose)\n",
    "Continuous_LinReg_features = ['Hometown_Distance','Age_Difference']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Create a dataset with our scaled data\n",
    "X_test_scaled = scaler.transform(merged_test_df[Continuous_LinReg_features])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Simply selecting the binary features of our data and converting into a numpy array\n",
    "X_test_binary = np.array(merged_test_df[Binary_LinReg_features])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Now we are merging our two datasets together into one\n",
    "X_test_penalized = np.concatenate((X_test_scaled, X_test_binary), axis = 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "PenalizedRegModel_Season25_Predictions = LassoCV_model.predict(X_test_penalized)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Ordinal Logistic Regression Predictions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "OrdLogRegModel_season25_Predictions = Ordinal_Logistic_Model.predict(merged_test_df[features])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Save and Submit your Predictions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Replace with the your model predictions where \"OrdLogRegModel_season25_Predictions\" currently is\n",
    "FINAL_PREDICTIONS = pd.DataFrame({'Name': merged_test_df['Contestant_Name'],\n",
    "                                  'ElimWeek': OrdLogRegModel_season25_Predictions})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "FINAL_PREDICTIONS"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You should have one prediction per contestant (we are excluding contestants who show up in the middle of the season from this contest). In total, your prediction dataset should have 32 rows and 2 columns:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "FINAL_PREDICTIONS.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Write the output to a csv file for your submission\n",
    "#PLEASE INCLUDE YOUR FIRST AND LAST NAME IN THE FILE NAME SO WE KNOW WHO TO GIVE PRIZES TO\n",
    "FINAL_PREDICTIONS.to_csv('Charlie_Putnam_predictions.csv', index = False)"
   ]
  }
 ],
 "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.7.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}