linear_regression.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "linear_regression.ipynb",
      "provenance": [],
      "collapsed_sections": [],
      "authorship_tag": "ABX9TyPI3u20XojxW7QsIFv0ctSN",
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/seanie12/linear-algebra/blob/master/linear_regression.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "L2KZJx68ZoAt",
        "colab_type": "text"
      },
      "source": [
        "# Linear Regression with Normal Equation and Gradient Descent"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "tNEN8ZdOZ1h-",
        "colab_type": "text"
      },
      "source": [
        "$X \\in \\mathfrak{M}_{m\\times n}(\\mathbb{R}), \\mathbf{y} \\in \\mathbb{R}^m, \\theta \\in \\mathbb{R}^n$ \n",
        "\\begin{equation}\n",
        "\\mathcal{L(\\theta) = || \\mathbf{X}\\theta - \\mathbf{y} ||^2}\n",
        "\\end{equation}"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "dcah4mjSZzVw",
        "colab_type": "text"
      },
      "source": [
        "We want to find $\\theta$ which minizes the error function $\\mathcal{L}(\\theta)$."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Yomto0rgZq2m",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "import numpy as np\n",
        "from numpy.linalg import inv"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "hH_yIDwMjCNK",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "num_features = 5\n",
        "X = np.random.randn(100, num_features)\n",
        "Y = np.random.randn(100)\n",
        "theta = np.random.randn(num_features)\n"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "OYZud5TsjC98",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 34
        },
        "outputId": "632d0699-77ba-4edf-bb90-56f0d07b2253"
      },
      "source": [
        "error = np.sum((Y - np.dot(X, theta))**2, axis=0)\n",
        "error"
      ],
      "execution_count": 3,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "441.1478116816043"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 3
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "P_ZyiamTjQXa",
        "colab_type": "text"
      },
      "source": [
        "## Normal Equation"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "GpOPQZ12jenJ",
        "colab_type": "text"
      },
      "source": [
        " It has closed-form solution as follows, which is Normal Equation :\n",
        "\\begin{align*}\n",
        "\\text{argmin}_{\\theta} \\mathcal{L}(\\theta) =  (X^tX)^{-1}X^t\\mathbf{y}\n",
        "\\end{align*}"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "RtjaJSRBjFlU",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 34
        },
        "outputId": "77783ba6-43a9-431a-9e2d-f3555fddc1f2"
      },
      "source": [
        "theta_hat = np.dot(inv(X.T @ X) @ X.T, Y)\n",
        "final_loss = np.sum((Y - np.dot(X, theta_hat))**2, axis=0)\n",
        "final_loss "
      ],
      "execution_count": 4,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "88.15259780536041"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 4
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ELMoTHFTjneC",
        "colab_type": "text"
      },
      "source": [
        "## Gradient Descent"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "yoipsqgjZNbM",
        "colab_type": "text"
      },
      "source": [
        "Since $\\mathcal{L}(\\theta)$ is differentiable function, we can find find the parameter $\\theta$ by gradient descent as follows, where $\\eta$ is learning rate:\n",
        "\\begin{align*}\n",
        "\\nabla_\\theta \\mathcal{L}(\\theta) := (\\frac{\\partial \\mathcal{L}(\\theta)}{\\partial \\theta_1}, \\ldots, \\frac{\\partial \\mathcal{L}(\\theta)}{\\partial \\theta_n}, )\n",
        "\\end{align*}\n",
        "\\begin{equation}\n",
        "\\theta^{(t+1)} = \\theta^{(t)} - \\eta \\nabla_{\\theta}\\mathcal{L}(\\theta^{(t)})\n",
        "\\end{equation}"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "6dLF4UNkfZOD",
        "colab_type": "text"
      },
      "source": [
        "Let $\\mathbf{x_i}$ the ith row vector of $X$.\n",
        "\\begin{align*}\n",
        "\\frac{\\partial \\mathcal{L}(\\theta)}{\\partial \\theta_j} &= \\sum\\limits_{i=1}^m (\\mathbf{x}_i\\cdot \\theta - y_i)^2 \\\\\n",
        "&= \\sum\\limits_{i=1}^m \\frac{\\partial (\\mathbf{x}_i \\cdot \\theta - y_i)^2}{\\partial \\theta_j} \\\\\n",
        "&= \\sum\\limits_{i=1}^m 2(\\mathbf{x}_i \\cdot \\theta - y_i) \\frac{\\partial (\\mathbf{x}_i \\cdot \\theta - y_i)}{\\partial \\theta_j} \\\\\n",
        "&=  \\sum\\limits_{i=1}^m 2(\\mathbf{x}_i \\cdot \\theta - y_i) X_{ij} \\:\\: \\text{ where } X_{ij} \\text{ is } (i,j) \\text{ entry of } X \\\\\n",
        "&= 2 [X]^j \\cdot (X\\theta - \\mathbf{y}) \\text{ where } [X]^j \\text{is the jth column vector of }X \\\\\n",
        "\\end{align*}\n",
        "\n",
        "\n",
        "\\begin{equation}\\therefore \\nabla_{\\theta}\\mathcal{L}(\\theta) = 2X^t(X\\theta - \\mathbf{y})\n",
        "\\end{equation}"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "rHLAhXYi8Fu5",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "class Regression(object):\n",
        "    def __init__(self, num_features, lr=0.001):\n",
        "        self.weight = np.random.randn(num_features)\n",
        "        self.lr = lr\n",
        "\n",
        "\n",
        "    def forward(self, x):\n",
        "        self.logits = np.dot(x, self.weight)\n",
        "        \n",
        "        return self.logits\n",
        "\n",
        "    def backward(self, x, y):\n",
        "        grad = 2 * x.T @ (self.logits - y)\n",
        "        self.weight = self.weight - self.lr * grad\n",
        "        return grad\n"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "3ONudBNb8F4H",
        "colab_type": "code",
        "outputId": "31c45a34-fe26-4d03-874c-3d2e94eff3ba",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 34
        }
      },
      "source": [
        "net = Regression(num_features, lr=0.001)\n",
        "\n",
        "for _ in range(100):\n",
        "    logits = net.forward(X)\n",
        "    loss = np.sum((logits - Y) ** 2, 0)\n",
        "    grad = net.backward(X, Y)\n",
        "\n",
        "logits = net.forward(X)\n",
        "sgd_loss = np.sum((logits - Y)**2, axis=0)\n",
        "print(sgd_loss)"
      ],
      "execution_count": 6,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "88.15259780539971\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "QljCgg0W8GAY",
        "colab_type": "code",
        "outputId": "f5196013-853f-4b5e-89e0-738fa33b0cfd",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 34
        }
      },
      "source": [
        "abs(final_loss - sgd_loss)"
      ],
      "execution_count": 7,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "3.929301328753354e-11"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 7
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "7JiYOz4P3k8Y",
        "colab_type": "text"
      },
      "source": [
        ""
      ]
    }
  ]
}