tutorials/multivariable_linear_regression_gradient_descent.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "fdb3af1b",
   "metadata": {},
   "source": [
    "## Gradient Descent and Multivariable Linear Regression"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8b7272e3",
   "metadata": {},
   "source": [
    "**Copyright (c) Meta Platforms, Inc. and affiliates.**\n",
    "\n",
    "This source code is licensed under the MIT license found in the\n",
    "LICENSE file in the root directory of this source tree."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c2514d65",
   "metadata": {},
   "source": [
    "In this notebook we are going to focus on performing a multivariable linear regression with Kotlin and DiffKt, specifically where we have two input variables and one output variable. \n",
    "\n",
    "A **linear regression** fits a linear function through data, often using a least squares method, to make predictions on new data. Using DiffKt we will perform gradient descent to optimize the coefficients. \n",
    "\n",
    "Bring in the DiffKt library for the tensor library and automatic differentiation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "d1b23168",
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "@file:DependsOn(\"../kotlin/api/build/libs/api.jar\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "be712fe4",
   "metadata": {},
   "source": [
    "We are going to need a few imports for this notebook. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "19db96c3",
   "metadata": {},
   "outputs": [],
   "source": [
    "import java.net.URL\n",
    "import org.diffkt.*\n",
    "import kotlin.random.Random"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3b53ebc9",
   "metadata": {},
   "source": [
    "Next up let's bring in our dataset. This dataset has two input variable columns `x1` and `x2` as well as an output variable column `y`. It only has 15 records. You can find the CSV data here: https://bit.ly/35ebET5.\n",
    "\n",
    "To start our model, use the `URL` and Kotlin `Sequence` to process the CSV into `Point` objects, which we also declare a class for. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "4072656f",
   "metadata": {},
   "outputs": [],
   "source": [
    "data class Point(val x1: Float, val x2: Float, val y: Float)\n",
    "\n",
    "val points = URL(\"https://bit.ly/35ebET5\")    // read CSV\n",
    "    .readText().split(Regex(\"\\\\r?\\\\n\"))       // split lines using regular expression\n",
    "    .filter { it.matches(Regex(\"[-,.0-9]+\")) }  // filter only numeric records using regular expression\n",
    "    .map { it.split(\",\").map{ it.toFloat()} } // split commas into columns\n",
    "    .map { (x1,x2,y) -> Point(x1,x2,y) }      // map to Point objects\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "becb2a02",
   "metadata": {},
   "source": [
    "Below is a visualization of the data in a 3D scatterplot. `x1` and `x2` are mapped to the horizontal axes, and the output variable `y` to the vertical axis. \n",
    "\n",
    "![](./resources/HiQXwPkaiO.mp4)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "faa563be",
   "metadata": {},
   "source": [
    "We are going to need to map these points to DiffKt tensors. We will use the `tensorOf()` function and map the `points` inside it. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "06c0ea25",
   "metadata": {},
   "source": [
    "Now when we map to the input tensor `x` and the output tensor `y`, we use lambda functions as arguments to specify what columns we want to generate and on what values. However notice on the `x` tensor below we add a third column simply returning a $ 1 $. This is going to add a column of 1's next to our `x1` and `x2` input variables. Why is this necessary? It will serve as a placeholder to generate the intercept coefficient. Without it, we would only generate the slopes for `x1` and `x2` without any intercept. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "2b9acb7d",
   "metadata": {},
   "outputs": [],
   "source": [
    "// map variables to input and output variable tensors\n",
    "// add a placeholder \"1\" column to generate intercept on input tensor\n",
    "val x = tensorOf(points.flatMap { listOf(it.x1, it.x2, 1f) }.map(::FloatScalar) ).reshape(points.size, 3)\n",
    "val y = tensorOf(points.map { it.y }.map(::FloatScalar) ).reshape(points.size, 1)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0262307a",
   "metadata": {},
   "source": [
    "Before we dive into the DiffKt operation and perform differentation and gradient descent, let's look at the mathematical side of things. \n",
    "\n",
    "A linear function with two inputs and one outputs will have three coefficients: $ \\beta_1 $ and $ \\beta_2 $ for the slopes of $ x_1 $ and $ x_2 $ respectively, and $ \\beta_0 $ for the intercept. Here is the linear function below. \n",
    "\n",
    "$ y = \\beta_2x_2 + \\beta_1x_1 + \\beta_0 $\n",
    "\n",
    "We need to solve for those three coefficients $ \\beta_0 $, $ \\beta_1 $, and $ \\beta_2 $. Rather than separate these as three separate scalar values separately, we can consolidate them into a single tensor holding three values. Let's initalize the `betas` tensors here. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "66152bcf",
   "metadata": {},
   "outputs": [],
   "source": [
    "// initialize coefficients\n",
    "var betas: DTensor = FloatTensor.random(Random,Shape(3,1))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6e0e4750",
   "metadata": {},
   "source": [
    "To visualize the tensor operations, let's say our betas were initialized with the following values. \n",
    "\n",
    "$ \\beta = \\left[\\begin{matrix}0.1\\\\0.2\\\\0.5\\end{matrix}\\right] $ \n",
    "\n",
    "And let's say we have 3 records of $ X $ inputs with the added column of 1's. \n",
    "\n",
    "$ X = \\left[\\begin{matrix}2 & 10 & 1\\\\4 & 20 & 1\\\\10 & 30 & 1\\end{matrix}\\right] $ \n",
    "\n",
    "To get the predicted `Y` values, we apply matrix multiplication (dot products) between the input $ X $ variables (with the additional column of 1's) and the $ \\beta $ coefficients. \n",
    "\n",
    "$ \\hat{Y} = X \\cdot \\beta $ \n",
    "\n",
    "$ \\hat{Y} = \\left[\\begin{matrix}2 & 10 & 1\\\\4 & 20 & 1\\\\10 & 30 & 1\\end{matrix}\\right] \\cdot \\left[\\begin{matrix}0.1\\\\0.2\\\\0.5\\end{matrix}\\right] $ \n",
    "\n",
    "$ \\hat{Y} =  \\left[\\begin{matrix}(2 \\times 0.1) + (10 \\times 0.2) + (1 \\times 0.5) \\\\(1 \\times 0.1) + (20 \\times 0.2) + (1 \\times 0.5) \\\\(10 \\times 0.1) + (30 \\times 0.2) + (1 \\times 0.5) \\end{matrix}\\right] $\n",
    "\n",
    "$ \\hat{Y} = \\left[\\begin{matrix}2.7\\\\4.9\\\\7.5\\end{matrix}\\right] $\n",
    "\n",
    "So that would yield predictions of $ 2.7 $, $ 4.9 $, and $ 7.5 $. \n",
    "\n",
    "To get predictions on all data given the current `betas` coefficients, use DiffKt's `*` operator: "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "8d86063b",
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[[0.068403006], [0.34552747], [-0.2825855], [-0.085838735], [1.207175], [0.82209784], [1.0211763], [0.52355623], [0.6804542], [0.09496629], [0.753953], [1.4067023], [0.38704538], [0.5864217], [0.40068573], [1.1661887], [0.09574813], [1.3083658], [1.0136129], [0.51993704], [0.0053132176], [0.4336069], [0.76465416], [0.20007932], [0.80502176], [0.12513003], [0.82067263], [0.36929575], [1.3200634], [0.078621864], [0.32192543], [0.22869426], [1.1282701], [0.0058833957], [0.7473586], [0.4677243], [0.7926725], [0.53712946], [0.40495625], [0.46880335], [0.58279735], [0.98099613], [1.2595377], [0.9098745], [0.07523024]]"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "val yPredictions = x.matmul(betas)\n",
    "yPredictions"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9fe81bf4",
   "metadata": {},
   "source": [
    "To calculate the total loss, let's use a sum of squared loss. Subtract the actual $ Y $ values from the predicted $ \\hat{Y} $ values. Take those differences, square them, and sum them. \n",
    "\n",
    "$ E = \\sum{(Y - \\hat{Y})^2 } $ \n",
    "\n",
    "Let's say we have these predicted $ \\hat{Y} $ and actual $ Y $ values. \n",
    "\n",
    "$ \\hat{Y} = \\left[\\begin{matrix}2.7\\\\4.9\\\\7.5\\end{matrix}\\right] $\n",
    "\n",
    "$ Y = \\left[\\begin{matrix}3.0\\\\5.0\\\\7.0\\end{matrix}\\right] $\n",
    "\n",
    "Here is how we would calculate the sum of squares. \n",
    "\n",
    "$ E = \\sum{(Y - \\hat{Y})^2 } $ \n",
    "\n",
    "$ E = \\sum{(\\left[\\begin{matrix}3.0\\\\5.0\\\\7.0\\end{matrix}\\right] - \\left[\\begin{matrix}2.7\\\\4.9\\\\7.5\\end{matrix}\\right])^2 } $ \n",
    "\n",
    "$ E = \\sum{(\\left[\\begin{matrix}0.3\\\\0.1\\\\-0.5\\end{matrix}\\right])^2} $ \n",
    "\n",
    "$ E = \\sum{\\left[\\begin{matrix}0.09\\\\0.01\\\\0.25\\end{matrix}\\right]} $ \n",
    "\n",
    "$ E = 0.35 $ \n",
    "\n",
    "We can implement this as a `loss()` function in Kotlin using DiffKt as shown below. Remember that the predicted $ \\hat{Y} $ values are the dot products of $ X $ and the $ \\beta $ coefficients. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "e2a13c8d",
   "metadata": {},
   "outputs": [],
   "source": [
    "fun loss(betas: DTensor): DScalar =\n",
    "    (y - (x.matmul(betas))).pow(2).sum()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "94b1dc45",
   "metadata": {},
   "source": [
    "Finally we are ready to perform gradient descent. For $ 100,000 $ iterations, we will use a learning rate of $ .001 $ and take the reverse derivative of the `loss()` function with regards to the `betas` tensor. This will return the gradient for each $ \\beta $ coefficient respectively which we multiply by the learning rate and subtract from the `betas` tensor. We subtract because we want to descend on the gradients."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "0455f74d",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "betas=[[0.3192737], [0.6086176], [-1.0212873]]"
     ]
    }
   ],
   "source": [
    "// The learning rate\n",
    "val lr = .001F\n",
    "\n",
    "// The number of iterations to perform gradient descent\n",
    "val iterations = 1_000\n",
    "\n",
    "// Perform gradient descent\n",
    "for (i in 0..iterations) {\n",
    "\n",
    "    // get gradients for line slope and intercept\n",
    "    val betaGradients = reverseDerivative(betas, ::loss)\n",
    "\n",
    "    // update m and b by subtracting the (learning rate) * (slope)\n",
    "    betas -= betaGradients * lr\n",
    "}\n",
    "print(\"betas=$betas\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2c1a0479",
   "metadata": {},
   "source": [
    "We should get the following coefficient values based on our `betas` tensors after this finishes training. \n",
    "\n",
    "$ \\beta_1 = 0.3192737 $ \n",
    "\n",
    "$ \\beta_2 = 0.6086176 $ \n",
    "\n",
    "$ \\beta_0 = -1.0212873 $ \n",
    "\n",
    "$ y = \\beta_1 x_1 +  \\beta_2 x_2  + \\beta_0 $ \n",
    "\n",
    "$ y = 0.3192737 x_1 +  0.6086176 x_2 - 1.0212873 $ \n",
    "\n",
    "If we were to visualize this linear regression fit as a plane in a 3D scatterplot, here is what it looks like. \n",
    "\n",
    "![](./resources/cCVILwOVBr.mp4)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Kotlin",
   "language": "kotlin",
   "name": "kotlin"
  },
  "language_info": {
   "codemirror_mode": "text/x-kotlin",
   "file_extension": ".kt",
   "mimetype": "text/x-kotlin",
   "name": "kotlin",
   "nbconvert_exporter": "",
   "pygments_lexer": "kotlin",
   "version": "1.7.0-dev-3303"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}