04b-time-series-prediction.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Sequential data modelling\n",
    "In this notebook we will work with the Jena weather dataset. \"The dataset is recorded at the Weather Station at the Max Planck Institute for Biogeochemistry in Jena, Germany.\" [source](https://www.kaggle.com/pankrzysiu/weather-archive-jena). The data is collected over the years 2009-2016. Using this dataset we will try to predict the temperature in the future.\n",
    "\n",
    "We will start by preprocessing the data, reshape it, scale it, and split it to train/val/test. The reshaping step is new and very common when dealing with sequential data.\n",
    "\n",
    "After preprocessing the data we will create a baseline model (dense network) and see how well it performs on this data.\n",
    "\n",
    "As usual we start by loading our data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "library(keras)\n",
    "library(ggplot2)\n",
    "source(\"04-helpers.R\")\n",
    "\n",
    "use_multi_cpu() # Remove this line if you are running in a local environment."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data <- load_jena_dataset() # This will download the dataset into a folder \"data\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Excercise 1\n",
    "As usual, explore the data a bit. What are the names of the columns/features? How many examples do we have?\n",
    "\n",
    "**Hint** these column names have spaces in the name so refering to them using `data$col a` will not work, instead do `data$\"col a\"`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We want to predict the temperature, `T (degC)`. Let us plot that feature and see how the temperature varies over those years."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ggplot(data, aes(x = 1:nrow(data),y=`T (degC)`)) + geom_line()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Excercise 2\n",
    "Looking at the plot above. How many years worth of data do we have? To answer this, just try counting the peaks. What do the peaks represent?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "YOUR ANSWER HERE"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "These are quite many measurements. In fact, a new measurement is drawn every 10 minutes, 6 measurements and hour, 6*24=144 measurements a day. Let us plot 10 days (144 * 10) of measurements. Starting from day 60.\n",
    "\n",
    "Feel free change the start_day or the duration."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "start_day = 60 \n",
    "duration_days = 10\n",
    "measurements_per_day = 6*24\n",
    "start_index = 1 + (start_day - 1) * measurements_per_day\n",
    "end_index = start_index + duration_days * measurements_per_day\n",
    "ggplot(data[start_index:end_index,], aes(x = start_index:end_index, y = `T (degC)`)) + geom_line()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Preprocessing"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Today we will need to spend some additional time preprocessing our data. This is because our sequential model will expect to get the data in dimensions `(examples, sequence_length, features)`. This will impact other preprocessing stages as well.\n",
    "\n",
    "Let us start with a basic step, removing features which we will not use and change our data object to a matrix."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Excercise 3.1\n",
    "Remove the `Date Time` feature from `data` and store the result in a variable called `features`. Verify that you have done this correctly by calling `names()` on `features` and check the `dim()`. The data should have dimension `(420551, 14)`.\n",
    "\n",
    "**Hint** `Date Time` is the first column in `data`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Excercise 3.2\n",
    "Now let us change `features` into a matrix. The type is important for our reshape step.. First check the type of `features` using the `typeof()` function. Then map `features` to a new variable `feature_matrix` by using the function `as.matrix()`. Verify that it has type `double`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Reshaping\n",
    "Now `feature_matrix` is of the correct type and contains all our features. Now we need to \"reshape\" this matrix (2D tensor) to an array (3D tensor). That is we need to change the dimensions from `(examples, features)` to `(examples, sequence_length, features)`. This is because sequence models expect the input to have shape `(examples, sequence_length, features)`. We need to create sequences out of multiple examples and we do this by reshaping our data. Let us see how the [`reshape()`](https://rstudio.github.io/reticulate/reference/array_reshape.html) function works."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Excercise 4.1\n",
    "Reshape the vector (1D tensor) `vec` to a `(2,2)` matrix (2D tensor)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "vec <- 1:4\n",
    "mat <- ?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The reshape function expects our inputs to fit perfectly into the new tensor. We can not reshape the vector above to a tensor with dimensions `(2,2,2)` as it requires 8 values."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Excercise 4.2\n",
    "Create a vector (1D tensor) with 8 elements (a sequence like above) and reshape it into a 3D tensor called `tensor`. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "vec <- ?\n",
    "tensor <- ?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Run the cell below to verify that your solution worked and notice the ordering of elements."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "tensor[1,,]\n",
    "tensor[2,,]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Going back to our task. Our data has the shape `(examples, features)` (2D tensor) and we need to reshape it to `(examples, sequence_length, features)` (3D tensor) by reducing the `examples` dimension. Run the example below to see how we might do this."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "example_1 <- 1:3\n",
    "example_2 <- 11:13\n",
    "example_3 <- 111:113\n",
    "example_4 <- 1111:1113\n",
    "x <- matrix(c(example_1, example_2, example_3, example_4), nrow = 4, byrow = TRUE)\n",
    "dim(x)\n",
    "x[1,]\n",
    "x[2,]\n",
    "x[3,]\n",
    "x[4,]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`x` represents our data but only has 4 examples, each with 3 features, `(4,3)`. Let us now reshape it to so that we we 2 examples, each with a sequence of length 2. That is, `(2,2,3)`. We do not change the feature dimension."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "reshaped_x <- array_reshape(x, c(2,2,3))\n",
    "reshaped_x[1,,]\n",
    "reshaped_x[2,,]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can still use the `apply` function on `reshaped_x`. The cell below will sum each feature over all the examples. This will help us scale our data later."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "apply(reshaped_x, 3, sum)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Success!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Excercise 4.3\n",
    "Reshape `feature_matrix` so that it becomes a 3D tensor with dimensions `(1752, 240, 14)`. By doing this we are create 1752 sequences of length 240. For this to work `feature_matrix` needs to fit into a 3D tensor, so only take the first 1752*240=420.480 elements from `feature_matrix`. Store the result in a variable called `feature_matrix_reshaped`. Verify the dimensions after reshaping."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sequence_length <- 240\n",
    "feature_matrix_reshaped <- ?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Wow! This was a lot of work, but well worth it, as understanding how the `array_reshape()` function works is very helpful."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Splitting\n",
    "Now we will split our `feature_matrix_reshaped` into train/val/test sets. We will, like usual, only split on the `examples` axis. To make the splits easy to understand and debug we will store 6 indices, one for where each section begins and ends."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Excercise 5.1\n",
    "Create the following 6 indices.\n",
    "- train_start_index\n",
    "- train_end_index\n",
    "- val_start_index\n",
    "- val_end_index\n",
    "- test_start_index\n",
    "- test_end_index\n",
    "\n",
    "Create these 6 indices so that:\n",
    "- `train_start_index = 1`\n",
    "- after the training set the validation set starts\n",
    "- after the validation set ends, the test set starts\n",
    "- `test_end_index = 1752`\n",
    "- Use the normal sizes of splits 0.8\\*0.8 / 0.8\\*0.2 / 0.2 for train/val/test.\n",
    "- Use `dim(feature_matrix_reshaped)[1]` to get the number of examples.\n",
    "\n",
    "**Hint** Be sure to return integers by using the function `as.integer( some_decimal )`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "train_start_index <- 1\n",
    "train_end_index <- ?\n",
    "val_start_index <- ?\n",
    "val_end_index <- ?\n",
    "test_start_index <- ?\n",
    "test_end_index <- ?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Run the cell below to verify that all the indices are correct, they should **roughly** be (it does not need to be exactly the same, off by 1 or 2 is alright):\n",
    "- train_start_index = 1\n",
    "- train_end_index = 1121\n",
    "- val_start_index = 1122\n",
    "- val_end_index = 1401\n",
    "- test_start_index = 1402\n",
    "- test_end_index = 1752\n",
    "\n",
    "If you are unable to do this, just store the indices as these numbers."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "train_start_index\n",
    "train_end_index\n",
    "val_start_index\n",
    "val_end_index\n",
    "test_start_index\n",
    "test_end_index"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Scaling\n",
    "Before splitting the `feature_matrix_reshaped` let us scale the features. Today we will scale all our features, even the feature we want to predict. This is because we will use previous values of that feature as input to our model, so we need to scale it. We will only use the training set to compute the mean and standard deviation and use that to scale the whole dataset. \n",
    "\n",
    "To only compute the mean and standard deviation on the training set we use the `apply()` function. Verify that this is done correctly in the cell below (not a trick question)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# We compute the mean and sd for each feature independently, The features are on the 3rd axis.\n",
    "mean <- apply(feature_matrix_reshaped[train_start_index:train_end_index,,], 3, mean)\n",
    "std <- apply(feature_matrix_reshaped[train_start_index:train_end_index,,], 3, sd)\n",
    "\n",
    "mean[2] # The mean temperature\n",
    "mean[1] # The mean air pressure"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This time we will not use the `scale()` function since it only works with (at most) 2D tensors, so we need to do this manually. We store the result in `feature_matrix_reshaped_scaled`. We provide the code for this, but we recommend reading through it."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create an empty 3D tensor with same dimensions as feature_matrix_reshaped.\n",
    "feature_matrix_reshaped_scaled <- array(0, dim = dim(feature_matrix_reshaped))\n",
    "\n",
    "# For each feature we subtract the mean and divide by the std.\n",
    "for (j in 1:dim(feature_matrix_reshaped)[3]) {\n",
    "    feature_matrix_reshaped_scaled[,,j] <- (feature_matrix_reshaped[,,j] - mean[j]) / std[j]\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let us see if this all makes sense and check the first temperature."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "feature_matrix_reshaped[1,1,2]\n",
    "feature_matrix_reshaped_scaled[1,1,2]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Brr... A cold day."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Excercise 5.2\n",
    "Now we can finish splitting the data into train/val/test sets. Create the 3 datasets `train_data`, `val_data` and `test_data` based on the indices above and the `feature_matrix_reshaped_scaled` data and verify the dimensions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "train_data <- ?\n",
    "val_data <- ?\n",
    "test_data <- ?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## The predictions\n",
    "Now we have successfully scaled and split our dataset. The last step is to create the correct predictions. For each sequences we need to create a corresponding target. To do this, we need to select a temperature point somewhere in the future and predict it. Let us try to predict the temperature for the next day, that is, 144 timesteps into the future."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To explain this further let us take an example. The first training sequence we have is `train_data[1,,]`. It contains 240 measurements, each made out of 14 features. We need to select the temperature 144 timesteps into the future from the last measurement in our sequence. This temperature is in the next training example, at `train_data[2,144,2]`.\n",
    "\n",
    "Let us now do this systematically for each dataset."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Excercise 6.1\n",
    "Create `y_train`, `y_val` and `y_test` from `train_data`, `val_data` and `test_data`, respectively. These will be our targets, `y`. Fill in the missing information to the cell below. Verify that the dimensions are correct (these will be vectors, so use the `length()` function). \n",
    "\n",
    "**Hint** There should be 1 fewer elements in the `y`'s since we cannot use the first sequence to use a prediction for anything. We will address this in the next excercise."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "y_train <- ?\n",
    "y_val <- ?\n",
    "y_test <- ?\n",
    "length(y_train)\n",
    "length(y_val)\n",
    "length(y_test)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since we needed to use the temperature from the next example to create our predictions, we would need to use the temperature from the first sequence in validation set as the correct prediction for the last sequence in the training set, and so on. We cannot be bothered to do that and we just skip the last sequence in the training data as input."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Excercise 6.2\n",
    "Create `x_train`, `x_val` and `x_test` from the `train_data`, `val_data` and `test_data`, respectively. Do not include the last sequence in `train_data` (and `val_data` and `test_data`) since we do not have the correct prediction for the last sequence in `y_train` (and `y_val` and `y_test`). Verify that `dim(x_train)[1]` is equal to `length(y_train)`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "x_train <- ?\n",
    "x_val <- ?\n",
    "x_test <- ?\n",
    "dim(x_train)[1]\n",
    "length(y_train)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We are finally done with all of the preprocessing. Wow! This was a lot of work! This is amount of preprocessing for sequential models is quite common as you want to be able to change many of the decisions we made in these steps. For example, we might want to increase/decrease timesteps in each sequence or predict longer/shorter into the future. This would not be possible if these preprocessing steps were done for us."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## The baseline model\n",
    "Now we want to create a dense network model to serve as our baseline. We then want to see the improvement we can get by adding an RNN.\n",
    "\n",
    "In this model notice a few things.\n",
    "- The dense networks expect input to be a vector, so the first layer in the network will flatten out the sequences (reversing a lot of our hard, reshape work).\n",
    "- Only a single hidden layer.\n",
    "- Since we are doing regression, the last layer has a single unit and no activation.\n",
    "- We use the mean absolute error, `mae`. It worked better."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model <- keras_model_sequential() %>%\n",
    "    layer_flatten(input_shape = c(sequence_length, 14)) %>%\n",
    "    layer_dense(units = 12, activation = 'relu') %>%\n",
    "    layer_dense(units = 1)\n",
    "\n",
    "cat(summary(model))\n",
    "\n",
    "model %>% compile(\n",
    "    optimizer = optimizer_adam(lr = 0.001),\n",
    "    loss = \"mae\"\n",
    ")\n",
    "\n",
    "history <- model %>% fit(\n",
    "    x = x_train,\n",
    "    y = y_train,\n",
    "    validation_data = list(x_val, y_val),\n",
    "    epochs = 80,\n",
    "    batch_size = 32,\n",
    "    callbacks=list(Progress$new())\n",
    ")\n",
    "\n",
    "plot(history)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Feel free to tweak the model a bit and then evaluate your model when you are done."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model %>% evaluate(x_test, y_test)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we do the same for a simple RNN network. In this network notice:\n",
    "- A single RRN layer.\n",
    "- Very few parameters, we use the same weights for every 240 timestep.\n",
    "- Fewer epochs.\n",
    "- Same loss.\n",
    "- Slow training (because of all the timesteps we have)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model <- keras_model_sequential() %>%\n",
    "    layer_simple_rnn(units = 12, input_shape = c(sequence_length, 14)) %>%\n",
    "    layer_dense(units = 1)\n",
    "\n",
    "cat(summary(model))\n",
    "\n",
    "model %>% compile(\n",
    "    optimizer = optimizer_adam(lr = 0.001),\n",
    "    loss = \"mae\"\n",
    ")\n",
    "\n",
    "history <- model %>% fit(\n",
    "    x = x_train,\n",
    "    y = y_train,\n",
    "    validation_data = list(x_val, y_val),\n",
    "    epochs = 40,\n",
    "    batch_size = 32,\n",
    "    callbacks=list(Progress$new())\n",
    ")\n",
    "\n",
    "plot(history)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Feel free to tweak the model a bit and then evaluate your model when you are done. You will notice that there is a lot more variation in the dense network than in the RNN."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model %>% evaluate(x_test, y_test)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we are finished! In this notebook worked with the Jena weather dataset. We first preprocessed the data using the reshape function, to create sequences. Then we scaled the data \"manually\" and finally split it to train/dev/test.\n",
    "\n",
    "We then established a baseline model, a dense network and then tried a simple RNN network. We saw that the RNN performed a lot more consistently, better, had way fewer parameters, but also took longer to train due to our long sequences."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Bonus excercise\n",
    "Scaling the temperature in the predictions, `y`, is not as nice as having them as the correct celsius values. Repeat the processing pipeline for `y` without scaling the predictions. Then run the same model using these predictions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Bonus excercise\n",
    "Try using the [`predict`](https://keras.rstudio.com/reference/predict.keras.engine.training.Model.html) function. This is the function you would run \"in production\" to get the actual predictions. Remember to shift and scale the output with the mean and std we computed earlier."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "R",
   "language": "R",
   "name": "ir"
  },
  "language_info": {
   "codemirror_mode": "r",
   "file_extension": ".r",
   "mimetype": "text/x-r-source",
   "name": "R",
   "pygments_lexer": "r",
   "version": "3.5.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}