notebooks/official/workbench/spark/spark_ml.ipynb

{
  "cells": [
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "ur8xi4C7S06n"
      },
      "outputs": [],
      "source": [
        "# Copyright 2022 Google LLC\n",
        "#\n",
        "# Licensed under the Apache License, Version 2.0 (the \"License\");\n",
        "# you may not use this file except in compliance with the License.\n",
        "# You may obtain a copy of the License at\n",
        "#\n",
        "#     https://www.apache.org/licenses/LICENSE-2.0\n",
        "#\n",
        "# Unless required by applicable law or agreed to in writing, software\n",
        "# distributed under the License is distributed on an \"AS IS\" BASIS,\n",
        "# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.\n",
        "# See the License for the specific language governing permissions and\n",
        "# limitations under the License."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "# SparkML with Dataproc and BigQuery\n",
        "\n",
        "<table align=\"left\">\n",
        "\n",
        "  <td>\n",
        "<a href=\"https://github.com/GoogleCloudPlatform/vertex-ai-samples/blob/main/notebooks/official/workbench/spark/spark_ml.ipynb\" target='_blank'>\n",
        "      <img src=\"https://cloud.google.com/ml-engine/images/github-logo-32px.png\" alt=\"GitHub logo\">\n",
        "      View on GitHub\n",
        "    </a>\n",
        "  </td>\n",
        "  <td>\n",
        "        <a href=\"https://colab.research.google.com/github/GoogleCloudPlatform/vertex-ai-samples/blob/main/notebooks/official/workbench/spark/spark_ml.ipynb\">\n",
        "        <img src=\"https://cloud.google.com/ml-engine/images/colab-logo-32px.png\" alt=\"Colab logo\"> Run in Colab\n",
        "        </a>\n",
        "  </td>\n",
        "  <td>\n",
        "<a href=\"https://console.cloud.google.com/vertex-ai/workbench/deploy-notebook?download_url=https://raw.githubusercontent.com/GoogleCloudPlatform/vertex-ai-samples/main/notebooks/official/workbench/spark/spark_ml.ipynb\" target='_blank'>\n",
        "      <img src=\"https://lh3.googleusercontent.com/UiNooY4LUgW_oTvpsNhPpQzsstV5W8F7rYgxgGBD85cWJoLmrOzhVs_ksK_vgx40SHs7jCqkTkCk=e14-rj-sc0xffffff-h130-w32\" alt=\"Vertex AI logo\">\n",
        "      Open in Vertex AI Workbench\n",
        "    </a>\n",
        "  </td>                                                                                               \n",
        "</table>"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "## Overview\n",
        "\n",
        "This notebook tutorial runs Apache SparkML jobs with Dataproc and BigQuery to exemplify a common machine learning pipeline use case: data ingestion and cleaning, feature engineering, modeling, and model evaluation.\n",
        "\n",
        "Learn more about [Vertex AI Workbench](https://cloud.google.com/vertex-ai/docs/workbench/introduction) and [Dataproc](https://cloud.google.com/vertex-ai/docs/pipelines/dataproc-component)."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "### Objective\n",
        "\n",
        "This tutorial runs an Apache SparkML job that fetches data from the BigQuery dataset, performs exploratory data analysis, cleans the data, executes feature engineering, trains the model, evaluates the model, outputs results, and saves the model to a Cloud Storage bucket.\n",
        "\n",
        "This tutorial uses the following Google Cloud ML services:\n",
        "\n",
        "- `Dataproc`\n",
        "- `BigQuery`\n",
        "- `Vertex AI Training`\n",
        "\n",
        "The steps performed are:\n",
        "\n",
        "- Sets up a Google Cloud project and Dataproc cluster.\n",
        "- Creates a Cloud Storage bucket and a BigQuery dataset.\n",
        "- Configures the spark-bigquery-connector.\n",
        "- Ingests BigQuery data into a Spark DataFrame.\n",
        "- Performa Exploratory Data Analysis (EDA).\n",
        "- Visualizes the data with samples.\n",
        "- Cleans the data.\n",
        "- Selects features.\n",
        "- Trains the model.\n",
        "- Outputs results.\n",
        "- Saves the model to a Cloud Storage bucket.\n",
        "- Deletes the resources created for the tutorial."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "### Dataset\n",
        "\n",
        "The [NYC TLC (Taxi and Limousine Commission) Trips](https://console.cloud.google.com/marketplace/product/city-of-new-york/nyc-tlc-trips) (New York taxi and limosine trips data) and [NYC Citi Bike Trips](https://console.cloud.google.com/marketplace/product/city-of-new-york/nyc-citi-bike) (NYC public bicycle sharing system data) datasets are available in [BigQuery Public Datasets](https://cloud.google.com/bigquery/public-data). BigQuery provides free querying of up to 1TB of data each month."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "### Costs \n",
        "\n",
        "This tutorial uses billable components of Google Cloud:\n",
        "\n",
        "* [Vertex AI](https://cloud.google.com/vertex-ai/pricing)\n",
        "* [Cloud Storage](https://cloud.google.com/storage/pricing)\n",
        "* [Dataproc](https://cloud.google.com/dataproc/pricing)\n",
        "* [BigQuery](https://cloud.google.com/bigquery/pricing)\n",
        "\n",
        "You can use the [Pricing Calculator](https://cloud.google.com/products/calculator/) to generate a cost estimate based on your projected usage."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "### Installation\n",
        "\n",
        "Install the following packages:"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "b5c0183fba17"
      },
      "outputs": [],
      "source": [
        "import os\n",
        "\n",
        "if os.getenv(\"IS_TESTING\"):\n",
        "    \"\"\"\n",
        "    Since the testing suite doesn't support testing on Dataproc clusters,\n",
        "    the testing environment is setup to replicate Dataproc via the following steps:\n",
        "    \"\"\"\n",
        "    JAVA_VER = \"8u332-b09\"\n",
        "    JAVA_FOLDER = \"/tmp/java\"\n",
        "    FILE_NAME = f\"openlogic-openjdk-{JAVA_VER}-linux-x64\"\n",
        "    TAR_FILE = f\"{JAVA_FOLDER}/{FILE_NAME}.tar.gz\"\n",
        "    DOWNLOAD_LINK = f\"https://builds.openlogic.com/downloadJDK/openlogic-openjdk/{JAVA_VER}/openlogic-openjdk-{JAVA_VER}-linux-x64.tar.gz\"\n",
        "    PYSPARK_VER = \"3.1.3\"\n",
        "\n",
        "    ! rm -rf $JAVA_FOLDER\n",
        "    ! mkdir $JAVA_FOLDER\n",
        "    # Download Open JDK 8. Spark requires Java to execute.\n",
        "    ! wget -P $JAVA_FOLDER $DOWNLOAD_LINK\n",
        "    os.environ[\"JAVA_HOME\"] = f\"{JAVA_FOLDER}/{FILE_NAME}\"\n",
        "    ! tar -zxf $TAR_FILE -C $JAVA_FOLDER\n",
        "    ! echo $JAVA_HOME\n",
        "\n",
        "    # Install PySpark and other packages.\n",
        "    ! pip install pyspark==$PYSPARK_VER geopandas pyarrow rtree seaborn numpy==1.19.5 -q"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "### Create a Dataproc cluster\n",
        "\n",
        "The Spark job executed in this notebook tutorial is compute intensive. Since the job can take a significant amount time to complete in a standard notebook environment, this notebook tutorial runs on a Dataproc cluster that is created with the Dataproc Component Gateway and Jupyter component installed on the cluster.\n",
        "\n",
        "**Existing Dataproc with Jupyter cluster?**: If you have a running Dataproc cluster that has the [Component Gateway and Jupyter component installed on the cluster](https://cloud.google.com/dataproc/docs/concepts/components/jupyter#gcloud-command), you can use it in this tutorial. If you plan to use it, skip this step, and go to `Switch your kernel`.\n",
        "\n",
        "Set and name and [compute region](https://cloud.google.com/compute/docs/regions-zones#available) for your new cluster. Your `CLUSTER_NAME` must be **unique within your Google Cloud project**. It must start with a lowercase letter, followed by up to 51 lowercase letters, numbers, and hyphens, and cannot end with a hyphen."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "4e3cf9e61673"
      },
      "outputs": [],
      "source": [
        "# The testing environment does not use a Dataproc cluster so cluster creation is skipped during testing.\n",
        "if not os.getenv(\"IS_TESTING\"):\n",
        "    CLUSTER_NAME = \"[your-cluster]\"  # @param {type: \"string\"}\n",
        "    CLUSTER_REGION = \"[your-region]\"  # @param {type: \"string\"}\n",
        "\n",
        "    if CLUSTER_REGION == \"[your-region]\":\n",
        "        CLUSTER_REGION = \"us-central1\"\n",
        "\n",
        "    print(f\"CLUSTER_NAME: {CLUSTER_NAME}\")\n",
        "    print(f\"CLUSTER_REGION: {CLUSTER_REGION}\")"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "ff30d97bf211"
      },
      "outputs": [],
      "source": [
        "if not os.getenv(\"IS_TESTING\"):\n",
        "    !gcloud dataproc clusters create $CLUSTER_NAME \\\n",
        "        --region=$CLUSTER_REGION \\\n",
        "        --enable-component-gateway \\\n",
        "        --image-version=2.0 \\\n",
        "        --optional-components=JUPYTER"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "#### Switch your kernel\n",
        "\n",
        "Your notebook kernel is listed at the top of the notebook page. Your notebook should run on the Python 3 kernel running on your Dataproc cluster.\n",
        "\n",
        "Select **Kernel > Change Kernel** from the top menu, then select `Python 3 on CLUSTER_NAME: Dataproc cluster in REGION (Remote)`."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "### Set your project ID\n",
        "\n",
        "**If you don't know your project ID**, try the following:\n",
        "* Run `gcloud config list`.\n",
        "* Run `gcloud projects list`.\n",
        "* See the support page: [Locate the project ID](https://support.google.com/googleapi/answer/7014113)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "set_project_id"
      },
      "outputs": [],
      "source": [
        "PROJECT_ID = \"[your-project-id]\"  # @param {type:\"string\"}\n",
        "\n",
        "# Set the project id\n",
        "! gcloud config set project {PROJECT_ID}"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "region"
      },
      "source": [
        "#### Region\n",
        "\n",
        "You can also change the `REGION` variable used by Vertex AI. Learn more about [Vertex AI regions](https://cloud.google.com/vertex-ai/docs/general/locations)."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "2dw8q9fdQEH5"
      },
      "outputs": [],
      "source": [
        "REGION = \"us-central1\""
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "7f1992d08421"
      },
      "outputs": [],
      "source": [
        "REGION = \"[your-region]\"  # @param {type: \"string\"}\n",
        "\n",
        "# If you do not specify a region, it is set to \"us-central1\".\n",
        "if not REGION or REGION == \"\" or REGION == \"[your-region]\":\n",
        "    REGION = \"us-central1\""
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "#### UUID\n",
        "\n",
        "To avoid name collisions you can create a UUID for the current notebook session, then append the UUID to the name of resources such as the Cloud Storage bucket and BigQuery dataset created in this tutorial."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "6988e7272299"
      },
      "outputs": [],
      "source": [
        "import random\n",
        "import string\n",
        "\n",
        "\n",
        "# Generate a uuid of a specifed length(default=8)\n",
        "def generate_uuid(length: int = 8) -> str:\n",
        "    return \"\".join(random.choices(string.ascii_lowercase + string.digits, k=length))\n",
        "\n",
        "\n",
        "UUID = generate_uuid()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "gcp_authenticate"
      },
      "source": [
        "### Authenticate your Google Cloud account\n",
        "\n",
        "Depending on your Jupyter environment, you may have to manually authenticate. Follow the relevant instructions below.\n",
        "\n",
        "**1. Vertex AI Workbench**\n",
        "* Do nothing as you are already authenticated.\n",
        "\n",
        "**2. Local JupyterLab instance, uncomment and run:**"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "ce6043da7b33"
      },
      "outputs": [],
      "source": [
        "# ! gcloud auth login"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "0367eac06a10"
      },
      "source": [
        "**3. Colab, uncomment and run:**"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "21ad4dbb4a61"
      },
      "outputs": [],
      "source": [
        "# from google.colab import auth\n",
        "# auth.authenticate_user()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "c13224697bfb"
      },
      "source": [
        "**4. Service account or other**\n",
        "* See how to grant Cloud Storage permissions to your service account at https://cloud.google.com/storage/docs/gsutil/commands/iam#ch-examples."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "bucket:mbsdk"
      },
      "source": [
        "### Create a Cloud Storage bucket\n",
        "\n",
        "Create a storage bucket to store intermediate artifacts such as datasets."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "bucket"
      },
      "outputs": [],
      "source": [
        "BUCKET_URI = f\"gs://your-bucket-name-{PROJECT_ID}-unique\"  # @param {type:\"string\"}"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "autoset_bucket"
      },
      "source": [
        "**Only if your bucket doesn't already exist**: Run the following cell to create your Cloud Storage bucket."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "91c46850b49b"
      },
      "outputs": [],
      "source": [
        "! gsutil mb -l {REGION} -p {PROJECT_ID} {BUCKET_URI}"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "## Tutorial\n",
        "\n",
        "### Import required libraries"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "pRUOFELefqf1"
      },
      "outputs": [],
      "source": [
        "import os\n",
        "\n",
        "import matplotlib.pyplot as plt\n",
        "import numpy as np\n",
        "import pandas as pd\n",
        "import seaborn as sns\n",
        "from geopandas import gpd\n",
        "from pyspark.ml.evaluation import RegressionEvaluator\n",
        "from pyspark.ml.feature import VectorAssembler\n",
        "from pyspark.ml.regression import GBTRegressor\n",
        "# A Spark Session is how you interact with Spark SQL to create Dataframes\n",
        "from pyspark.sql import SparkSession\n",
        "# PySpark functions\n",
        "from pyspark.sql.functions import (col, floor, pandas_udf, to_timestamp, udf,\n",
        "                                   unix_timestamp)\n",
        "# These allow us to create a schema for our data\n",
        "from pyspark.sql.types import BooleanType, DoubleType\n",
        "from shapely.geometry import Point"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "71ecd729c4b9"
      },
      "source": [
        "**Note**: After importing libraries, if you see `ERROR 1: PROJ: proj_create_from_database: Open of /opt/conda/miniconda3/share/proj failed`, you can ignore it. This is due to a bug with a `geopandas` dependency.\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "init_aip:mbsdk,all"
      },
      "source": [
        "### Initialize the SparkSession\n",
        "\n",
        "To use Apache Spark with BigQuery, you must include the [spark-bigquery-connector](https://github.com/GoogleCloudDataproc/spark-bigquery-connector) when you initialize the SparkSession."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "init_aip:mbsdk,all"
      },
      "outputs": [],
      "source": [
        "VER = \"0.26.0\"\n",
        "FILE_NAME = f\"spark-bigquery-with-dependencies_2.12-{VER}.jar\"\n",
        "\n",
        "if os.getenv(\"IS_TESTING\"):\n",
        "    # The local testing environment is not configured to read GCS files using Spark. This is a workaround for testing.\n",
        "    connector = f\"https://github.com/GoogleCloudDataproc/spark-bigquery-connector/releases/download/{VER}/{FILE_NAME}\"\n",
        "else:\n",
        "    connector = f\"gs://spark-lib/bigquery/{FILE_NAME}\"\n",
        "\n",
        "\n",
        "# Initialize the SparkSession.\n",
        "spark = (\n",
        "    SparkSession.builder.appName(\"spark-ml-taxi-citibike\")\n",
        "    .config(\"spark.jars\", connector)\n",
        "    .getOrCreate()\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "init_aip:mbsdk,all"
      },
      "source": [
        "### Fetch data from BigQuery\n",
        "\n",
        "This tutorial uses the 2017 NYC Taxi dataset and the NYC Citi Bike dataset from 2013 to 2018. Since the time range of each dataset differs, match both datasets to the 2017 time period."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "init_aip:mbsdk,all"
      },
      "outputs": [],
      "source": [
        "# Load NYC_taxi in Github Activity Public Dataset from BigQuery.\n",
        "taxi_df = (\n",
        "    spark.read.format(\"bigquery\")\n",
        "    .option(\"table\", \"bigquery-public-data.new_york_taxi_trips.tlc_green_trips_2017\")\n",
        "    .load()\n",
        ")\n",
        "\n",
        "\n",
        "# Load NYC_Citibike in Github Activity Public dataset from BQ.\n",
        "bike_df = (\n",
        "    spark.read.format(\"bigquery\")\n",
        "    .option(\"table\", \"bigquery-public-data.new_york_citibike.citibike_trips\")\n",
        "    .load()\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "init_aip:mbsdk,all"
      },
      "source": [
        "Use the `to_timestamp()` and `unix_timestamp()` functions to convert `starttime` and `stoptime` to a Unix timestamp type.\n",
        "\n",
        "A Unix timestamp, which is an integer that represents time as the number of seconds since `January 1st, 1970 at 00:00:00 UTC`, is easy to parse and use across different systems."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "0156fd0e2921"
      },
      "outputs": [],
      "source": [
        "# Represents 2017-01-01 00:00:00\n",
        "START_2017 = 1483228800\n",
        "# Represents 2017-12-31 23:59:59\n",
        "END_2017 = 1514764799\n",
        "\n",
        "# Convert the type of starttime from a string to a Unix timestamp.\n",
        "bike_df = bike_df.withColumn(\n",
        "    \"starttime\", unix_timestamp(to_timestamp(col(\"starttime\")))\n",
        ")\n",
        "\n",
        "# Convert the type of stoptime from a string to a Unix timestamp.\n",
        "bike_df = bike_df.withColumn(\"stoptime\", unix_timestamp(to_timestamp(col(\"stoptime\"))))\n",
        "\n",
        "# Filter data so that `bike_df` only contains 2017 trips.\n",
        "bike_df = bike_df.where(\n",
        "    (col(\"starttime\") >= START_2017)\n",
        "    & (col(\"starttime\") <= END_2017)\n",
        "    & (col(\"stoptime\") >= START_2017)\n",
        "    & (col(\"stoptime\") <= END_2017)\n",
        ")\n",
        "\n",
        "\"\"\"\n",
        "Since notebook runtime scales linearly with the size of data, to limit processing time and compute resource consumption, \n",
        "this tutorial uses 20% of the Citi Bike dataset. You can try different values or remove the next line to use all of the data. \n",
        "\"\"\"\n",
        "bike_df = bike_df.sample(0.2)\n",
        "\n",
        "if os.getenv(\"IS_TESTING\"):\n",
        "    taxi_df = taxi_df.sample(0.0001)\n",
        "    bike_df = bike_df.sample(0.0001)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "### Perform Exploratory Data Analysis (EDA)\n",
        "\n",
        "As a first step, use Exploratory Data Analysis (EDA) to find patterns, relations, and anomalies in the data using graphs and other visualizations.\n",
        "\n",
        "Check Taxi dataset data types first."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "adeaf84f4054"
      },
      "outputs": [],
      "source": [
        "taxi_df.printSchema()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "Filter out unnecessary columns, and check null counts of the fields. Note that `pickup_location_id` and `dropoff_location_id` were converted from the BigQuery dataset into string format. To use it effectively, reformat these values."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "e5072e605966"
      },
      "outputs": [],
      "source": [
        "taxi_df = taxi_df.select(\n",
        "    \"pickup_datetime\",\n",
        "    \"dropoff_datetime\",\n",
        "    \"trip_distance\",\n",
        "    \"fare_amount\",\n",
        "    col(\"pickup_location_id\").cast(\"int\").alias(\"start_zone_id\"),\n",
        "    col(\"dropoff_location_id\").cast(\"int\").alias(\"end_zone_id\"),\n",
        ")\n",
        "\n",
        "taxi_df.describe().show()\n",
        "taxi_df.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "From this summary, you are able to know a lot of information.\n",
        "\n",
        "  * There are over 11 million trip histories for Green Taxi in 2017, which include some abnormal values, such as negative values.\n",
        "  * Until 2016, exact latitude and longitude were used for pickup and the dropoff locations. This data raised [privacy concerns](https://agkn.wordpress.com/2014/09/15/riding-with-the-stars-passenger-privacy-in-the-nyc-taxicab-dataset/). `pickup_location_id` and `dropoff_location_id` in the dataset corresponded to the [NYC Taxi Zones](https://data.cityofnewyork.us/Transportation/NYC-Taxi-Zones/d3c5-ddgc), roughly based on NYC Department of City Planning’s Neighborhood Tabulation Areas (NTAs), and are meant to indicate approximate neighborhood locations.\n",
        "\n",
        "Based on numbers shown above, filter out unrealistic values.\n",
        "\n",
        "  * Remove null values.\n",
        "  * trip_distance must be in between 0 and 20 (miles) since Manhattan is approximately 13.4 miles long.\n",
        "  * fare_amount, which does not include taxes and tips, must be in between 0 and 500 (dollars)."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "38aa94651f34"
      },
      "outputs": [],
      "source": [
        "taxi_df = taxi_df.where(\n",
        "    (col(\"trip_distance\") > 0)\n",
        "    & (col(\"trip_distance\") < 20)\n",
        "    & (col(\"fare_amount\") > 0)\n",
        "    & (col(\"fare_amount\") < 500)\n",
        ").dropna()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "Use the `unix_timestamp()` function to convert `pickup_datetime` and `dropoff_datetime` as a Unix Timestamp type.\n",
        "\n",
        "After converting `pickup_datetime` and `dropoff_datetime` to `start_time` and `end_time`, you can manipulate these timestamps produce more insteresting columnnar data."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "8be30eddd676"
      },
      "outputs": [],
      "source": [
        "# Convert the type of datetime from a string to a Unix timestamp.\n",
        "taxi_df = taxi_df.withColumn(\"start_time\", unix_timestamp(col(\"pickup_datetime\")))\n",
        "taxi_df = taxi_df.withColumn(\"end_time\", unix_timestamp(col(\"dropoff_datetime\")))\n",
        "\n",
        "# Calculate if this is a weekday\n",
        "taxi_df = taxi_df.withColumn(\n",
        "    \"is_weekday\",\n",
        "    ((floor(col(\"start_time\") / 86400) + 4) % 7 > 0)\n",
        "    & ((floor(col(\"start_time\") / 86400) + 4) % 7 < 6),\n",
        ")\n",
        "\n",
        "# Convert start_time to start_time_in_minute.\n",
        "taxi_df = taxi_df.withColumn(\n",
        "    \"start_time_in_hour\", floor((col(\"start_time\") % 86400) / 60 / 60)\n",
        ")\n",
        "\n",
        "# Calculate trip_duration.\n",
        "taxi_df = taxi_df.withColumn(\"trip_duration\", col(\"end_time\") - col(\"start_time\"))"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "### Perform Feature Engineering\n",
        "\n",
        "The Taxi dataset contains trips for all NYC boroughs, while the Citi Bike dataset has trips between Manhattan and parts of Brooklyn and the Bronx. To match locations, filter out all boroughs except Manhattan.\n",
        "\n",
        "To convert the location ID to the borough, download the GeoJSON format of NYC Taxi zones from NYC OpenData."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "b94d6eaf7870"
      },
      "outputs": [],
      "source": [
        "# Load the GeoJSON format of NYC Taxi zones. If a URLError is thrown, please re-run this cell.\n",
        "gdf_zone = gpd.read_file(\n",
        "    \"https://data.cityofnewyork.us/api/geospatial/d3c5-ddgc?method=export&format=GeoJSON\"\n",
        ")"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "2be027cd9619"
      },
      "outputs": [],
      "source": [
        "# Convert the type of \"location_id\" column.\n",
        "gdf_zone[\"location_id\"] = gdf_zone[\"location_id\"].astype(\"long\")\n",
        "\n",
        "# Convert GeoPandas format to a dictionary.\n",
        "hm = gdf_zone.to_dict(\"index\")\n",
        "\n",
        "# Create a set that contains location_id if the location is in Manhattan.\n",
        "manhattan_set = {hm[i][\"location_id\"] for i in hm if hm[i][\"borough\"] == \"Manhattan\"}\n",
        "\n",
        "\n",
        "@udf(returnType=BooleanType())\n",
        "def is_in_manhattan(location_id):\n",
        "    \"\"\"\n",
        "    The preprocessing function takes location_id and returns whether the given location_id is in the manhattan_set.\n",
        "    Args:\n",
        "        location_id: A number from 0 to 263 that represents NYC Taxi Zone.\n",
        "    Returns:\n",
        "        A Boolean value of whether the given location_id is in the manhattan_set.\n",
        "    \"\"\"\n",
        "    return location_id in manhattan_set"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "10215f2acf7b"
      },
      "outputs": [],
      "source": [
        "# Trips occur within Manhattan and between different taxi zones.\n",
        "taxi_df = taxi_df.where(\n",
        "    (col(\"start_zone_id\") != col(\"end_zone_id\"))\n",
        "    & (is_in_manhattan(col(\"start_zone_id\")))\n",
        "    & (is_in_manhattan(col(\"end_zone_id\")))\n",
        ")"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "70511564352c"
      },
      "outputs": [],
      "source": [
        "taxi_df.show()\n",
        "taxi_df.describe().show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "Prior to visualizing this modified Taxi dataset, do similar feature engineering for the Citi Bike dataset."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "1f7b3436a287"
      },
      "outputs": [],
      "source": [
        "bike_df.printSchema()\n",
        "bike_df.describe().show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "This dataset reveals interesting information:\n",
        "\n",
        "  * There are over 2.8 million logs Citi Bike trips in 2017.\n",
        "  * The current dataset has some abnormal values.\n",
        "  * `starttime` and `stoptime` are in string format. To use them effectively, reformat them.\n",
        "  * The Citi Bike dataset does not contain the `location_id` for NYC Taxi Zones; only coordinates data. To convert the coordinates to the NYC Taxi Zones, use Pandas UDF, a Vectorized UDF feature (see [Introducing Pandas UDF for PySpark](https://www.databricks.com/blog/2017/10/30/introducing-vectorized-udfs-for-pyspark.html))."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "2be027cd9619"
      },
      "outputs": [],
      "source": [
        "@pandas_udf(\"long\")\n",
        "def preprocess_zone_id(lat: pd.Series, lon: pd.Series) -> pd.Series:\n",
        "    \"\"\"\n",
        "    The preprocessing function takes latitudes and longitudes and returns the corresponding Taxi Zone id.\n",
        "    Args:\n",
        "        lat: The latitude of the position.\n",
        "        lon: The longitude of the position.\n",
        "    Returns:\n",
        "        The Taxi Zone id.\n",
        "    \"\"\"\n",
        "    point_var = [Point(xy) for xy in zip(lon, lat)]\n",
        "    gdf_points = gpd.GeoDataFrame(\n",
        "        pd.DataFrame({\"lat\": lat, \"lon\": lon}), crs=\"epsg:4326\", geometry=point_var\n",
        "    )\n",
        "    gdf_joined = gpd.sjoin(gdf_points, gdf_zone, how=\"left\")\n",
        "    return gdf_joined[\"location_id\"]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "77e00e95bf23"
      },
      "outputs": [],
      "source": [
        "# Rename \"tripduration\" to \"trip_duration\".\n",
        "bike_df = bike_df.withColumnRenamed(\"tripduration\", \"trip_duration\")\n",
        "\n",
        "# Check whether the starttime is a weekday or a weekend.\n",
        "bike_df = bike_df.withColumn(\n",
        "    \"is_weekday\",\n",
        "    ((floor(col(\"starttime\") / 86400) + 4) % 7 > 0)\n",
        "    & ((floor(col(\"starttime\") / 86400) + 4) % 7 < 6),\n",
        ")\n",
        "\n",
        "# Convert starttime to start_time_in_minute\n",
        "bike_df = bike_df.withColumn(\n",
        "    \"start_time_in_hour\", floor((col(\"starttime\") % 86400) / 60 / 60)\n",
        ")\n",
        "\n",
        "# Add \"start_zone_id\" and \"end_zone_id\" by applying the UDF \"preprocess_zone_id\".\n",
        "bike_df = bike_df.withColumn(\n",
        "    \"start_zone_id\",\n",
        "    preprocess_zone_id(\n",
        "        bike_df[\"start_station_latitude\"], bike_df[\"start_station_longitude\"]\n",
        "    ),\n",
        ")\n",
        "bike_df = bike_df.withColumn(\n",
        "    \"end_zone_id\",\n",
        "    preprocess_zone_id(\n",
        "        bike_df[\"end_station_latitude\"], bike_df[\"end_station_longitude\"]\n",
        "    ),\n",
        ")\n",
        "\n",
        "# Filter `bike_df` using the following conditions:\n",
        "# 1) Trips must done between different taxi zones.\n",
        "# 2) The start time must be earlier than the end time.\n",
        "# 3) Trips must start and end within Manhattan.\n",
        "bike_df = bike_df.where(\n",
        "    (col(\"start_zone_id\") != col(\"end_zone_id\"))\n",
        "    & (col(\"starttime\") < col(\"stoptime\"))\n",
        "    & (is_in_manhattan(col(\"start_zone_id\")))\n",
        "    & (is_in_manhattan(col(\"end_zone_id\")))\n",
        ")\n",
        "\n",
        "bike_df = bike_df.select(\n",
        "    \"trip_duration\",\n",
        "    \"starttime\",\n",
        "    \"stoptime\",\n",
        "    \"start_station_latitude\",\n",
        "    \"start_station_longitude\",\n",
        "    \"end_station_latitude\",\n",
        "    \"end_station_longitude\",\n",
        "    \"usertype\",\n",
        "    \"start_zone_id\",\n",
        "    \"end_zone_id\",\n",
        "    \"is_weekday\",\n",
        "    \"start_time_in_hour\",\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "Assuming that Citi Bike users move through Manhattan's streets, which are the perpendicular, the trip distance can be calculated by applying [Manhattan distance](https://en.wikipedia.org/wiki/Taxicab_geometry).\n",
        "\n",
        "Manhattan distance calculation example: Assume you are at **W 15th St and 9th Ave, where Google NYC is located,** and you are going to **W 33rd St and 5th Ave, where the Empire State building is located**. You move east from 15th Street until you reach 5th Avenue, and then you travel north until you reach W 33rd St. In this case, **W 16th St and 5th Ave** or **W 33rd St and 9th Ave** can be a hinge point. If we set the starting point as **S**, the ending point as **E**, and the hinge point as **H**, the formula of Manhattan distance is `distance(S, H) + distance(H, E)`.\n",
        "\n",
        "However, additional calculation is needed to obtain the \"real\" Manhattan distance. First, unlike the Cartesian coordinate system, actual distance in longitude varies with latitude. For example, 1 degree of longitude represents roughly 69 miles at the equator, while it represents approximately 49 miles at 45 degrees north or south latitude. To calculate the exact distance in the real world using latitudes and longitudes, use [Haversine distance](https://en.wikipedia.org/wiki/Haversine_formula).\n",
        "\n",
        "Also, Manhattan streets are inclined at about 29 degrees from true north, so points should be rotateds S and E by 29 degrees counter-clockwise.\n",
        "\n",
        "Therefore, the trip distance formula should be `manhattan_distance = haversine_distance(S', H') + haversine_distance(H', E')` where `S'`, `H'`, `E'` are the rotated points."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "82b6e83a1d26"
      },
      "outputs": [],
      "source": [
        "@udf(returnType=DoubleType())\n",
        "def manhattan_dist(lat_start, lon_start, lat_end, lon_end):\n",
        "    \"\"\"\n",
        "    The preprocessing function takes the latitude and longitude of start position and end position,\n",
        "    and returns the Manhattan distance between them.\n",
        "    Args:\n",
        "        lat_start: The latitude of the start position.\n",
        "        lon_start: The longitude of the start position.\n",
        "        lat_end: The latitude of the end position.\n",
        "        lon_end: The longitude of the end position.\n",
        "    Returns:\n",
        "        The Manhattan distance between start and end position.\n",
        "    \"\"\"\n",
        "    # Approximate radius of Earth in mile.\n",
        "    EARTH_RADIUS = 3958.76\n",
        "\n",
        "    # Approximate inclined degree of the streets in Manhattan.\n",
        "    THETA = np.radians(-28.904)\n",
        "\n",
        "    def haversine(pos_start, pos_end):\n",
        "        \"\"\"\n",
        "        The helper function takes start and end position and returns the haversine distance.\n",
        "        Args:\n",
        "            pos_start: a latitude and a longitude of start position in np.array form.\n",
        "            pos_end: a latitude and a longitude of end position in np.array form.\n",
        "        Returns:\n",
        "            The Haversine, the spherical distance between `pos_start` and `pos_end` on a sphere, in this case, the earth.\n",
        "        \"\"\"\n",
        "        rad_dist_lat = np.radians(pos_end[0] - pos_start[0])\n",
        "        rad_dist_lon = np.radians(pos_end[1] - pos_start[1])\n",
        "        rad_start_lat = np.radians(pos_start[0])\n",
        "        rad_end_lat = np.radians(pos_end[0])\n",
        "        dist = 2 * np.arcsin(\n",
        "            np.sqrt(\n",
        "                np.sin(rad_dist_lat / 2) ** 2\n",
        "                + np.cos(rad_start_lat)\n",
        "                * np.cos(rad_end_lat)\n",
        "                * np.sin(rad_dist_lon / 2) ** 2\n",
        "            )\n",
        "        )\n",
        "        return EARTH_RADIUS * dist\n",
        "\n",
        "    def rotate(pos, theta):\n",
        "        \"\"\"\n",
        "        The helper function takes position and the degree theta, and returns rotated position.\n",
        "        Args:\n",
        "            pos: a latitude and a longitude of position in np.array form.\n",
        "            theta: the degree to rotate.\n",
        "        Returns:\n",
        "            Rotated position.\n",
        "        \"\"\"\n",
        "        rotate = np.array(\n",
        "            [[np.cos(theta), np.sin(theta)], [-np.sin(theta), np.cos(theta)]]\n",
        "        )\n",
        "        return np.matmul(rotate, pos)\n",
        "\n",
        "    if not (lat_start and lon_start and lat_end and lon_end):\n",
        "        return -1\n",
        "\n",
        "    # Convert positions to np.array format.\n",
        "    start, end = np.array([lat_start, lon_start]), np.array([lat_end, lon_end])\n",
        "\n",
        "    # Rotate each positions by 29' using a helper function.\n",
        "    rotated_start, rotated_end = rotate(start, THETA), rotate(end, THETA)\n",
        "\n",
        "    # Get rotated hinge point using rotated start and end point.\n",
        "    rotated_hinge = np.array([rotated_start[0], rotated_end[1]])\n",
        "\n",
        "    # Re-rotate the hinge point.\n",
        "    hinge = rotate(rotated_hinge, -THETA)\n",
        "\n",
        "    # Return the sum of the Haversine distance between start and hinge, hinge and end.\n",
        "    return float(haversine(start, hinge) + haversine(hinge, end))"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "7ece73e72bea"
      },
      "outputs": [],
      "source": [
        "# Calculate the Manhattan distance between \"start_station\" and \"end_station\".\n",
        "bike_df = bike_df.withColumn(\n",
        "    \"trip_distance\",\n",
        "    manhattan_dist(\n",
        "        col(\"start_station_latitude\"),\n",
        "        col(\"start_station_longitude\"),\n",
        "        col(\"end_station_latitude\"),\n",
        "        col(\"end_station_longitude\"),\n",
        "    ),\n",
        ")"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "bc5dd05670e6"
      },
      "outputs": [],
      "source": [
        "# Filter out the bike_df with the same condition as the taxi_df\n",
        "bike_df = bike_df.where(\n",
        "    (col(\"trip_distance\") > 0)\n",
        "    & (col(\"trip_distance\") < 20)\n",
        "    & (col(\"trip_duration\") > 0)\n",
        "    & (col(\"start_station_latitude\") > 40)\n",
        "    & (col(\"start_station_latitude\") < 41)\n",
        "    & (col(\"end_station_latitude\") > 40)\n",
        "    & (col(\"end_station_latitude\") < 41)\n",
        "    & (col(\"start_station_longitude\") > -75)\n",
        "    & (col(\"start_station_longitude\") < -73)\n",
        "    & (col(\"end_station_longitude\") > -75)\n",
        "    & (col(\"end_station_longitude\") < -73)\n",
        ")"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "29d2a801cf6a"
      },
      "outputs": [],
      "source": [
        "bike_df.printSchema()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "Check the distributions for the numerical columns. In PySpark, visualizing is expensive since the data is large. \n",
        "\n",
        "For example, the NYC Taxi dataset has more than 10M rows. Therefore, approximately 3% of total data (approximately 330k rows) are extracted as a sample."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "42a319b0c66c"
      },
      "outputs": [],
      "source": [
        "if not os.getenv(\"IS_TESTING\"):\n",
        "    taxi_pd = taxi_df.sample(0.03).toPandas()\n",
        "else:\n",
        "    taxi_pd = taxi_df.toPandas()\n",
        "\n",
        "# Convert \"object\" type columns to the float type.\n",
        "taxi_pd[\"trip_distance\"] = taxi_pd[\"trip_distance\"].astype(float)\n",
        "taxi_pd[\"fare_amount\"] = taxi_pd[\"fare_amount\"].astype(float)\n",
        "taxi_pd.info()\n",
        "\n",
        "if not os.getenv(\"IS_TESTING\"):\n",
        "    bike_pd = bike_df.sample(0.03).toPandas()\n",
        "else:\n",
        "    bike_pd = bike_df.toPandas()\n",
        "bike_pd.info()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "Remove extreme data points, and bin the data into a box and histogram."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "723d2f8fa827"
      },
      "outputs": [],
      "source": [
        "TAXI_COLUMNS_TO_SHOW = {\"trip_distance\", \"trip_duration\", \"fare_amount\"}\n",
        "\n",
        "taxi_pd_filtered = taxi_pd.query(\n",
        "    \"trip_distance > 0 and trip_distance < 20 \\\n",
        "    and trip_duration > 0 and trip_duration < 10000\"\n",
        ")\n",
        "\n",
        "sns.relplot(\n",
        "    data=taxi_pd_filtered,\n",
        "    x=\"trip_distance\",\n",
        "    y=\"trip_duration\",\n",
        "    kind=\"scatter\",\n",
        ")\n",
        "\n",
        "for column in taxi_pd_filtered.columns:\n",
        "    if column in TAXI_COLUMNS_TO_SHOW:\n",
        "        _, ax = plt.subplots(1, 2, figsize=(10, 4))\n",
        "        taxi_pd_filtered[column].plot(kind=\"box\", ax=ax[0])\n",
        "        taxi_pd_filtered[column].plot(kind=\"hist\", ax=ax[1])\n",
        "        plt.title(column)\n",
        "        plt.figure()\n",
        "plt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "Most columns are skewed right, and the `trip_distance` increases as `trip_duration` increases. Most trips are completed within an hour (3600 sec).\n",
        "\n",
        "Do a similar analysis for the Citi Bike dataset."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "355716305c5f"
      },
      "outputs": [],
      "source": [
        "BIKE_COLUMNS_TO_SHOW = {\"trip_distance\", \"trip_duration\"}\n",
        "\n",
        "bike_pd_filtered = bike_pd.query(\n",
        "    \"trip_distance > 0 and trip_distance < 20 \\\n",
        "    and trip_duration > 0 and trip_duration < 10000\"\n",
        ")\n",
        "\n",
        "sns.relplot(\n",
        "    data=bike_pd_filtered,\n",
        "    x=\"trip_distance\",\n",
        "    y=\"trip_duration\",\n",
        "    hue=\"usertype\",\n",
        "    kind=\"scatter\",\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "The plot for `trip_distance` and `trip_duration` for Citi Bike is different than the same plot for the Taxi dataset.\n",
        "\n",
        "Many trip durations are longer than expected given the trip distance. The apparent reasons for this behavior are 1) unlike taxis riders, city bike riders are less likely to have a specific trip destination, and 2) occasional Citi Bike users don't park their bicycles at the station. \n",
        "\n",
        "To remove outliers, let's assume the average bike speed is at least 2 miles per hour. Let's filter out trips where the `trip_duration` multiplied by speed is less than `trip_distance`, and let's remove the `Customer` user type."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "45112e073997"
      },
      "outputs": [],
      "source": [
        "bike_pd_filtered = bike_pd.query(\n",
        "    \"trip_distance > 0 and trip_distance < 8 \\\n",
        "    and (trip_duration * 0.0005) <= trip_distance \\\n",
        "    and trip_duration > 0 and trip_duration < 10000 \\\n",
        "    and usertype=='Subscriber'\"\n",
        ")\n",
        "\n",
        "sns.relplot(\n",
        "    data=bike_pd_filtered,\n",
        "    x=\"trip_distance\",\n",
        "    y=\"trip_duration\",\n",
        "    hue=\"usertype\",\n",
        "    kind=\"scatter\",\n",
        ")\n",
        "\n",
        "for column in bike_pd_filtered.columns:\n",
        "    if column in BIKE_COLUMNS_TO_SHOW:\n",
        "        _, ax = plt.subplots(1, 2, figsize=(10, 4))\n",
        "        bike_pd_filtered[column].plot(kind=\"box\", ax=ax[0])\n",
        "        bike_pd_filtered[column].plot(kind=\"hist\", ax=ax[1])\n",
        "        plt.title(column)\n",
        "        plt.figure()\n",
        "plt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "Based on the visualized data, apply a filter to remove unrealistic values."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "09d3908a6012"
      },
      "outputs": [],
      "source": [
        "taxi_df = taxi_df.where(\n",
        "    (col(\"trip_duration\") > 0)\n",
        "    & (col(\"trip_duration\") < 3600)\n",
        "    & ((col(\"trip_duration\") * 0.00055) <= col(\"trip_distance\"))\n",
        ")\n",
        "\n",
        "bike_df = bike_df.where(\n",
        "    (col(\"trip_duration\") > 0)\n",
        "    & (col(\"trip_duration\") < 3600)\n",
        "    & (col(\"usertype\") == \"Subscriber\")\n",
        "    & ((col(\"trip_duration\") * 0.00055) <= col(\"trip_distance\"))\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "### Feature Selection\n",
        "\n",
        "Not all features in our dataset will be useful. \n",
        "Since the purpose of this tutorial is to compare two datasets, we need to use the same features for training. There are some useful features for training the Taxi dataset, such as `fare_amount`, but since the Citi Bike dataset does not have it, it will not be used in this tutorial.\n",
        "\n",
        "After choose the following columns as features, they must be assembled with `VectorAssembler()`, which is feature transformer that merges multiple columns into a vector column."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "1ca926fe3b73"
      },
      "outputs": [],
      "source": [
        "feature_cols = [\n",
        "    \"is_weekday\",\n",
        "    \"start_time_in_hour\",\n",
        "    \"start_zone_id\",\n",
        "    \"end_zone_id\",\n",
        "    \"trip_distance\",\n",
        "]\n",
        "\n",
        "assembler = VectorAssembler(inputCols=feature_cols, outputCol=\"features\")\n",
        "\n",
        "# Transform each column into the vector form.\n",
        "taxi_transformed_data = assembler.transform(taxi_df)\n",
        "bike_transformed_data = assembler.transform(bike_df)\n",
        "\n",
        "# Split randomly with training and test data.\n",
        "(taxi_training_data, taxi_test_data) = taxi_transformed_data.randomSplit([0.95, 0.05])\n",
        "(bike_training_data, bike_test_data) = bike_transformed_data.randomSplit([0.95, 0.05])"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "### Training the Model\n",
        "\n",
        "PySpark provides several Regression Models, including  `LinearRegression`, `GeneralizedLinearRegression`, `DecisionTreeRegressor`, `RandomForestRegressor`, and `GBTRegressor`.\n",
        "`GBTRegressor`, known as `Gradient Boosted Trees` is a popular regression method using ensembles of decision trees. This algorithm iteratively trains decision trees in order to minimize a loss function.\n",
        "Although this model is computationally expensive, testing showed that `GBT Regressor` provided the best results."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "4d8ef0874cbf"
      },
      "outputs": [],
      "source": [
        "# Define GBTRegressor\n",
        "gbt = GBTRegressor(\n",
        "    featuresCol=\"features\",\n",
        "    labelCol=\"trip_duration\",\n",
        "    predictionCol=\"pred_trip_duration\",\n",
        ")\n",
        "\n",
        "# Define evaluator for r2 score.\n",
        "evaluator_r2 = RegressionEvaluator(\n",
        "    labelCol=gbt.getLabelCol(), predictionCol=gbt.getPredictionCol(), metricName=\"r2\"\n",
        ")\n",
        "\n",
        "# Define evaluator for RMSE error.\n",
        "evaluator_rmse = RegressionEvaluator(\n",
        "    labelCol=gbt.getLabelCol(), predictionCol=gbt.getPredictionCol(), metricName=\"rmse\"\n",
        ")"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "c0ecfdae3db2"
      },
      "outputs": [],
      "source": [
        "# Train a model and get predictions for the Taxi dataset. It takes several minutes.\n",
        "taxi_gbt_model = gbt.fit(taxi_training_data)\n",
        "taxi_gbt_predictions = taxi_gbt_model.transform(taxi_test_data)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "e43b652e6107"
      },
      "outputs": [],
      "source": [
        "# Train a model and get predictions for the Citi Bike dataset. It takes several minutes.\n",
        "bike_gbt_model = gbt.fit(bike_training_data)\n",
        "bike_gbt_predictions = bike_gbt_model.transform(bike_test_data)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "d4197f357a7b"
      },
      "outputs": [],
      "source": [
        "# Evaluate the r2 score for the Taxi dataset\n",
        "taxi_gbt_accuracy_r2 = evaluator_r2.evaluate(taxi_gbt_predictions)\n",
        "print(f\"Taxi Test GBT R2 Accuracy = {taxi_gbt_accuracy_r2}\")\n",
        "\n",
        "# Evaluate the RMSE for the Taxi dataset\n",
        "taxi_gbt_accuracy_rmse = evaluator_rmse.evaluate(taxi_gbt_predictions)\n",
        "print(f\"Taxi Test GBT RMSE Accuracy = {taxi_gbt_accuracy_rmse}\")"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "410b54ce8541"
      },
      "outputs": [],
      "source": [
        "# Evaluate the r2 score for the Citi Bike dataset\n",
        "bike_gbt_accuracy_r2 = evaluator_r2.evaluate(bike_gbt_predictions)\n",
        "print(f\"Bike Test GBT R2 Accuracy = {bike_gbt_accuracy_r2}\")\n",
        "\n",
        "# Evaluate the RMSE for the Citi Bike dataset\n",
        "bike_gbt_accuracy_rmse = evaluator_rmse.evaluate(bike_gbt_predictions)\n",
        "print(f\"Bike Test GBT RMSE Accuracy = {bike_gbt_accuracy_rmse}\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "### View the result\n",
        "\n",
        "The trained model for the Taxi dataset shows about 77% for r2 score and 230 root-mean square error (RMSE), while the model for the Citi Bike dataset shows about 71% for r2 score and 250 RMSE. There are some reasons for these relatively lower scores.\n",
        "\n",
        "* Didn't use the CrossValidator\n",
        "\n",
        "  * In real machine learning projects, cross-validation is an essential tool that allows us to better utilize data. However, it reruns the training algorithm several times, so it requres a signficant amount of time and computation resources to complete. For more information, see [Cross-validation (statistics)](https://en.wikipedia.org/wiki/Cross-validation_%28statistics%29).\n",
        "\n",
        "* Didn't use each feature in the Taxi dataset.\n",
        "\n",
        "  * The Taxi dataset has a `fare_amount` column, which has a strong correlation between the label, `trip_duration`. If that column is included, the model for the Taxi dataset shows 90+% r2 score.\n",
        "  * The reason why Taxi dataset had limited features is that the purpose of this project is to compare Taxi and the Citi Bike.\n",
        "  \n",
        "* Didn't use the exact coordinates for both datasets.\n",
        "\n",
        "  * Due to privacy concerns, the Taxi dataset no longer provides the precise locations of pickup and dropoff locations, only NYC Taxi zone locations.\n",
        "  * Taxi Zone is a categorical number that does not represents a geographical zone defined by NYC Taxi."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XoEqT2Y4DJmf"
      },
      "source": [
        "### Save the model to a Cloud Storage path"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "a2280850da83"
      },
      "outputs": [],
      "source": [
        "BUCKET_FOLDER = \"/tmp/bucket\"\n",
        "\n",
        "# In the testing environment, saving to Cloud Storage bucket routes to the local file system.\n",
        "if os.getenv(\"IS_TESTING\"):\n",
        "    ! rm -rf $BUCKET_FOLDER\n",
        "    ! mkdir $BUCKET_FOLDER\n",
        "    bike_gbt_model.write().overwrite().save(f\"{BUCKET_FOLDER}\")\n",
        "    ! gsutil cp -r $BUCKET_FOLDER $BUCKET_URI\n",
        "else:\n",
        "    bike_gbt_model.write().overwrite().save(f\"{BUCKET_URI}/\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "TpV-iwP9qw9c"
      },
      "source": [
        "## Cleaning up\n",
        "\n",
        "See [Clean up](https://cloud.google.com/vertex-ai/docs/workbench/managed/create-managed-notebooks-instance-console-quickstart#clean-up) to delete your project or the managed notebook created in this tutorial."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "TpV-iwP9qw9c"
      },
      "source": [
        "### Delete Google Cloud Storage bucket"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "d7a1a1a5e978"
      },
      "outputs": [],
      "source": [
        "import os\n",
        "\n",
        "DELETE_BUCKET = True\n",
        "\n",
        "if DELETE_BUCKET or os.getenv(\"IS_TESTING\"):\n",
        "    ! gsutil rm -r $BUCKET_URI"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "TpV-iwP9qw9c"
      },
      "source": [
        "### Delete Dataproc Cluster\n",
        "\n",
        "It is not possible to delete the cluster you are currently using unless you switch the kernel to local. To delete it, you need to switch the kernel to local `Python 3` or `PySpark`, set your `CLUSTER_NAME` and `CLUSTER_REGION` manually in the following cell, and execute the `gcloud` command.\n",
        "\n",
        "See [Deleting a cluster](https://cloud.google.com/dataproc/docs/guides/manage-cluster#console) to delete the Dataproc cluster created in this tutorial."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "5f6cae539ccc"
      },
      "outputs": [],
      "source": [
        "CLUSTER_NAME = \"[your-cluster-name]\"\n",
        "CLUSTER_REGION = \"[your-cluster-region]\""
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "61244eb74c17"
      },
      "outputs": [],
      "source": [
        "import os\n",
        "\n",
        "if not os.getenv(\"IS_TESTING\"):\n",
        "    ! gcloud dataproc clusters delete $CLUSTER_NAME --region=$CLUSTER_REGION -q"
      ]
    }
  ],
  "metadata": {
    "colab": {
      "collapsed_sections": [],
      "name": "spark_ml.ipynb",
      "toc_visible": true
    },
    "kernelspec": {
      "display_name": "Python 3",
      "name": "python3"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}