04. Polynomial Regression - R.ipynb

{
  "cells": [
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "Exercise 4 - Polynomial Regression\n===\n\nSometimes our data doesn't have a linear relationship, but we still want to predict an outcome.\n\nSuppose we want to predict how satisfied people might be with a piece of fruit. We would expect satisfaction would be low if the fruit is under-ripe or over-ripe, and satisfaction would be high if fruit is ripe, i.e. between under-ripe and over-ripe.\n\nThis is not something linear regression will help us with, however polynomial regression can help us make predictions for more complex non-linear relationships such as these.\n\nStep 1\n---\n\nIn this exercise, we will look at a dataset analysing internet traffic over the course of the day. Observations were made every hour, on the hour, over the course of several days. Suppose we want to predict the level of internet traffic we might see at any time during the day, how might we do this?\n\nLet's start by loading the required libraries for this session, and loading our data to have a look at it.\n\n** In the cell below replace: **\n\n** 1. `<dataPreviewStr>` with `str` **\n\n** 2. `<dataPreviewHead>` with `head` **\n\n** then __run the code__. **"
    },
    {
      "metadata": {
        "trusted": false
      },
      "cell_type": "code",
      "source": "# Load required libraries\nsuppressMessages(library(\"tidyverse\"))\n\ntraf_by_hr <- read.csv(\"Data/traffic_by_hour.csv\", check.names = FALSE)\n\n###\n# REPLACE <dataPreviewStr> WITH str AND <dataPreviewHead> WITH head\n###\n<dataPreviewStr>(traf_by_hr)\n<dataPreviewHead>(traf_by_hr)\n###",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "By observing the structure of `traf_by_hr`, we can observe that we have the following variables: \n\n* `Hour` 00 - 23, spread across the columns;\n* `Observation` 1 - 6, with each observation representing one row;\n* `Traffic` (in units Gbps), representing the values of each observation made every hour. \n\nStep 2\n---\n\nNext, we need to reshape the data using the package `tidyr` so we can plot our data using `ggplot2`. The `ggplot2` functions expect our data input to be in \"long\" format, i.e. a column for each feature name, and a row for each observation. Currently, our data is stored in \"wide\" format, with one feature, i.e. the hour, spread across multiple columns. We need to reshape the `traf_by_hr` data so that our variables, `Hour`, `Observation`, and `Traffic` are the column names.\n\n> \"Long\" and \"wide\" format are visual metaphors that describe two ways of presenting the exact same information. Hopefully you will understand the metaphor once you see \"long\" and \"wide\" data for yourself, as per below!\n\n#### Run the following code block to reshape `traf_by_hr` (\"wide\" format) to `traf_by_hr_long` (\"long\" format). You do not need to edit the code block."
    },
    {
      "metadata": {
        "trusted": false
      },
      "cell_type": "code",
      "source": "# Run this!\n\n# DO NOT EDIT\n\n# Reshape data to long format using tidyr's gather function\ntraf_by_hr_tall <- traf_by_hr %>%\ngather(key = \"Hour\", value = \"Traffic\") %>%\nmutate(Observation = as.factor(rep(1:6, 24)))  %>% \nmutate(Hour = as.numeric(Hour))\n\n# Check structure of reshaped data\nstr(traf_by_hr_tall)\nhead(traf_by_hr_tall, n = 10)",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "Compare the structural difference between `traf_by_hr_long` and `traf_by_hr` for yourself, particularly the dimensions. You should observe the following:\n\n* `traf_by_hr` (\"wide\" data): 6 observations (rows) x 24 variables (columns);\n* `traf_by_hr_long` (\"long\" data): 144 observations (rows) x 3 variables (columns).\n\nNow with our data in \"long\" format, we can use `ggplot2` functions for plotting.\n\nIn the `ggplot` function call within the code block below:\n\n** In the cell below replace: **  \n** 1. `<xData>` with `Hour` **  \n** 2. `<yData>` with `Traffic` **  \n** then __run the code__. **"
    },
    {
      "metadata": {
        "trusted": false
      },
      "cell_type": "code",
      "source": "# Visualise the data using ggplot2\ntraf_by_hr_tall  %>% \n\n###\n# REPLACE THE <xData> TO Hour AND <yData> TO Traffic\n###\nggplot(aes(x = <xData>, y = <yData>, colour = Observation)) +\n###\ngeom_line() +\nggtitle(\"Internet traffic for each hour of the day\") +\nxlab(\"Hour of the day (24 hour time)\") +\nylab(\"Internet traffic (Gbps)\") +\ntheme(plot.title = element_text(hjust = 0.5))",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "This plot looks a bit busy due to overplotting. We should summarize the data to help us visualize trends.\n\nStep 3\n---\n\nLet's see if we can visualize a clearer pattern by taking the __average values__ for each hour.\n\nIn the code block below, find the `geom_point` function call, then:\n\n** In the cell below replace:** \n\n**1. `<xData>` with `Hour`**  \n\n**2. `<yData>` with `Traffic_mean`**  \n\n**then __run the code__. **"
    },
    {
      "metadata": {
        "trusted": false
      },
      "cell_type": "code",
      "source": "# Find mean values for each hour using `dplyr` function `mutate`\ntraf_by_hr_tall <- traf_by_hr_tall  %>% \ngroup_by(Hour)  %>% \nmutate(Traffic_mean = mean(Traffic))  %>% \nungroup()  %>% \nas.data.frame()\n\n# Check structure of data\nstr(traf_by_hr_tall)\nhead(traf_by_hr_tall, n = 10)\n\n# Create plot\ntraf_by_hr_tall  %>% \nggplot() +\n\n# Plot the average of the 6 observations as points\n\n###\n# REPLACE <xData> TO Hour and <yData> to Traffic_mean\n###\ngeom_point(aes(x = <xData>, y = <yData>)) +\n###\n# Plot each observation 1 - 6 as a line\ngeom_line(aes(x = Hour, y = Traffic, colour = Observation), alpha = 0.5) +\nylab(\"Average internet traffic (Gbps)\") +\ntheme(plot.title = element_text(hjust = 0.5))",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "The plot above shows the average value for each hour as points (black), together with observations 1 - 6 as lines.\n\nWe can also plot our data using a graph type that summarizes the data for us, such as a box and whisker plot.\n\nIn the code below, within the `ggplot` function call, change the x and y variables to the features we want to observe.\n\n** In the cell below replace:**  \n**1. `<xData>` with `Hour`**  \n**2. `<yData>` with `Traffic`**  \n**then __run the code__. **"
    },
    {
      "metadata": {
        "trusted": false
      },
      "cell_type": "code",
      "source": "# Create box and whisker plot\ntraf_by_hr_tall  %>% \n###\n# REPLACE <xData> TO Hour and <yData> to Traffic\n###\nggplot(aes(x = <xData>, y = <yData>, group = Hour)) +\n###\ngeom_boxplot() +\nggtitle(\"Internet traffic for each hour of the day\") +\nxlab(\"Hour of the day (24 hour time)\") +\nylab(\"Internet traffic (Gbps)\") +\n# Align title to centre\ntheme(plot.title = element_text(hjust = 0.5))",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "This summarization of the data could help us make a prediction if we wanted to know the expected traffic exactly on the hour.\n\nBut, we'll need to be a bit clever if we want to make a good prediction of times in between."
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "Step 4\n---\n\nLet's use the midpoints in between the hours to help us analyse the relationship between the time of day and the amount of internet traffic.\n\nThe `lm` (linear model) function together with the `poly` (polynomial) function allow us to do just this. We need to specify a feature $x$ (time of day), our outcome $y$ (the amount of internet traffic), and the $degree$ of the polynomial (how curvy the line is).\n\n> You can use the `lm` function directly, but for this exercise we will use `lm` indirectly through `ggplot2`.\n\nFirst we will test polynomial functions with degrees 1, 2, 3 and 4. Note the first degree polynomial (degree = 1) has already been completed for you; first degree polynomials are linear, so we will include it for comparison.\n\n** In the cell below replace: **  \n** 1. `<addDegree2>` with `2`**  \n** 2. `<addDegree3>` with `3`**  \n** 3. `<addDegree4>` with `4`**  \n** then __run the code__.**"
    },
    {
      "metadata": {
        "trusted": false
      },
      "cell_type": "code",
      "source": "traf_by_hr_tall  %>% \nggplot(aes(x = Hour, y = Traffic_mean)) +\ngeom_point(alpha = 0.5) +\nstat_smooth(method = \"lm\", formula = y ~ poly(x, degree = 1), colour = \"black\", linetype = \"dashed\", se = FALSE) +\n\n###\n# REPLACE <addDegree2> WITH 2, <addDegree3> WITH 3, AND <addDegree4> WITH 4\n###\nstat_smooth(method = \"lm\", formula = y ~ poly(x, degree = <addDegree2>), colour = \"#F8766D\", se = FALSE) + # red\nstat_smooth(method = \"lm\", formula = y ~ poly(x, degree = <addDegree3>), colour = \"#00BFC4\", se = FALSE) + # blue\nstat_smooth(method = \"lm\", formula = y ~ poly(x, degree = <addDegree4>), colour = \"#7CAE00\", se = FALSE) + # green\n\nxlab(\"Hour of the day (24 hour time)\") +\nylab(\"Internet traffic (Gbps)\") +\nggtitle(\"Testing fit of polynomial functions (degrees 1 - 4) to internet traffic data\") +\ntheme(plot.title = element_text(hjust = 0.5))",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "None of these polynomial functions do a great job of generalising the data. Let's try a few more.\n\nIn the code below, test polynomial functions of degree 5, 6 and 7.\n\n** In the cell below replace:**  \n** 1. `<addDegree5>` with `5`**  \n** 2. `<addDegree6>` with `6`**  \n** 3. `<addDegree7>` with `7`**  \n** then __run the code__.**"
    },
    {
      "metadata": {
        "trusted": false
      },
      "cell_type": "code",
      "source": "traf_by_hr_tall  %>% \nggplot(aes(x = Hour, y = Traffic_mean)) +\ngeom_point(alpha = 0.5) +\n# Polynomials of degree 1 are linear, keep this for comparison\nstat_smooth(method = \"lm\", formula = y ~ poly(x, degree = 1), colour = \"black\", linetype = \"dashed\", se = FALSE) +\n# Change degree = ? to degree = 5\n\n###\n# REPLACE <addDegree5> WITH 5, <addDegree6> WITH 6, AND <addDegree7> WITH 7\n###\nstat_smooth(method = \"lm\", formula = y ~ poly(x, degree = <addDegree5>), colour = \"#F8766D\", se = FALSE) + # red\nstat_smooth(method = \"lm\", formula = y ~ poly(x, degree = <addDegree6>), colour = \"#00BFC4\", se = FALSE) + # blue\nstat_smooth(method = \"lm\", formula = y ~ poly(x, degree = <addDegree7>), colour = \"#7CAE00\", se = FALSE) + # green\n###\n\nxlab(\"Hour of the day (24 hour time)\") +\nylab(\"Internet traffic (Gbps)\") +\nggtitle(\"Testing fit of polynomial functions (degrees 5 - 7) to internet traffic data\") +\ntheme(plot.title = element_text(hjust = 0.5))",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "It looks like the 6th and 7th degree polynomials have an identical curve, so either of these polynomials will be a good model to use.\n\nWe could use an even higher degree polynomial to fit the model to our data even more tightly, but we don't want to overfit the curve, since we just want a generalization of the relationship between time of day and internet traffic.\n\nLet's see how our 6th degree polynomial alone compares to the real data.\n\n#### Run the code below."
    },
    {
      "metadata": {
        "trusted": false
      },
      "cell_type": "code",
      "source": "# Run this!\n\n# DO NOT EDIT\n\ntraf_by_hr_tall  %>% \nggplot(aes(x = Hour, y = Traffic, colour = Observation)) +\ngeom_line(alpha = 0.5) +\nstat_smooth(method = \"lm\", formula = y ~ poly(x, degree = 6), colour = \"black\", se = FALSE) +\nxlab(\"Hour of the day (24 hour time)\") +\nylab(\"Internet traffic (Gbps)\") +\nggtitle(\"Testing fit of 6th degree polynomial to internet traffic data\") +\ntheme(plot.title = element_text(hjust = 0.5))",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "Looking good!\n\nStep 5\n---\n\nNow let's try using this polynomial regression model to make a prediction for how much internet traffic there will be at a certain time of day. Let's choose the time 12:30pm.\n\n** Replace the `<addTime>` with `12.5`, which represents the time 12.30 pm, and run the code. **"
    },
    {
      "metadata": {
        "trusted": false
      },
      "cell_type": "code",
      "source": "###\n# REPLACE THE <addTime> WITH 12.5 (THIS REPRESENTS 12:30PM)\n###\nt <- data.frame(Hour = <addTime>)\n###\n\n# Save our 6th degree polynomial model\nlm_poly_6th <- lm(formula = Traffic ~ poly(Hour, degree = 6), data = traf_by_hr_tall)\n\n# Use predict function, and round the result to 2 decimal places\n# Input to predict function must be data frame, with column name set as x value\nt_pred <- round(predict(lm_poly_6th, t), 2)\n\nprint(paste(\"Based on our polynomial regression model, at time t = 12.5\",\n      \"the expected internet traffic is\", t_pred , \"Gbps.\"))\n\ntraf_by_hr_tall %>% \nggplot(aes(x = Hour, y = Traffic)) +\nstat_smooth(method = \"lm\", formula = y ~ poly(x, degree = 6), colour = \"black\", se = FALSE) +\n# Show predicted value as a point in red\ngeom_point(x = 12.5, y = t_pred, size = 3) +\n# Add horizontal line\ngeom_hline(yintercept = t_pred, linetype = \"dashed\", colour = \"red\") +\ngeom_vline(xintercept = 12.5, linetype = \"dashed\", colour = \"red\") +\nggtitle(\"Prediction of expected traffic at t = 12.5\") +\nxlab(\"Hour of the day (24 hour time)\") +\nylab(\"Internet traffic (Gbps)\") +\ntheme(plot.title = element_text(hjust = 0.5))",
      "execution_count": null,
      "outputs": []
    },
    {
      "metadata": {},
      "cell_type": "markdown",
      "source": "Conclusion\n---\n\nThere we have it! You have made a polynomial regression model and used it for analysis! This models gives us a prediction for the level of internet traffic we should expect to see at any given time of the day.\n\nYou can go back to the course and either click __'Next Step'__ to start an optional step with tips on how to better work with AI models, or you can go to the next module where instead of predicting numbers we predict categories."
    }
  ],
  "metadata": {
    "kernelspec": {
      "name": "r",
      "display_name": "R",
      "language": "R"
    },
    "language_info": {
      "mimetype": "text/x-r-source",
      "name": "R",
      "pygments_lexer": "r",
      "version": "3.4.1",
      "file_extension": ".r",
      "codemirror_mode": "r"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 2
}