Vectorizing_Gradient_Descent_Multivariate_Linear_Regression.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this article, I shall go over the topic of arriving at the **Vectorized Gradient-Descent formulae for the Cost function of the for Matrix form of training-data Equations.** And along with that the Fundamentals of Calculus (especially Partial Derivative) and Matrix Derivatives necessary to understand the process.\n",
    "\n",
    "**So our target of this article is to understand the full Mathematics and the flow behind arriving at the below formulae**, which is the Vectorized Gradient of the training-data Matrix\n",
    "\n",
    "![](https://imgur.com/wutr4wJ.png)\n",
    "\n",
    "\n",
    "### First a Refresher on basic Matrix Algebra\n",
    "\n",
    "A matrix A over a field K or, simply, a matrix A (when K is implicit) is a rectangular array of scalars usually presented in the following form:\n",
    "\n",
    "![](https://imgur.com/nprQ9nC.png)\n",
    "\n",
    "\n",
    "The rows of such a matrix A are the m horizontal lists of scalars:\n",
    "\n",
    "![](https://imgur.com/lm1dvCf.png)\n",
    "\n",
    "\n",
    "and the columns of A are the n vertical lists of scalars:\n",
    "\n",
    "![](https://imgur.com/PURIYkm.png)\n",
    "\n",
    "\n",
    "A matrix with m rows and n columns is called an m by n matrix, written m*n. The pair of numbers m and n is called the size of the matrix. Two matrices A and B are equal, written A = B, if they have the same size and if corresponding elements are equal. Thus, the equality of two m * n matrices is equivalent to a system of mn equalities, one for each corresponding pair of elements.\n",
    "A matrix with only one row is called a row matrix or row vector, and a matrix with only one column is called a column matrix or column vector. A matrix whose entries are all zero is called a zero matrix and will usually be denoted by 0.\n",
    "\n",
    "### Matrix Multiplication\n",
    "\n",
    "The below image is taken from Khan Academy’s excellent linear algebra course.\n",
    "\n",
    "![](https://imgur.com/D8o4Vau.png)\n",
    "\n",
    "**In above, each entry in the product matrix is the dot product of a row in the first matrix and a column in the second matrix**\n",
    "\n",
    "**More explanation for higher dimension case** — If the product **AB = C** is defined, where C is denoted by [cij], then the\n",
    "element cij is obtained by multiplying the elements in the ith row of **A** by the corresponding elements in the jth column of **B** and adding. Thus, if **A** has order k * n, and **B** has order n * p then\n",
    "\n",
    "![](https://imgur.com/6UbO5Ga.png)\n",
    "\n",
    "\n",
    "**There are four simple rules that will help us in multiplying matrices, listed here**\n",
    "\n",
    "1. Firstly, we can only multiply two matrices when the number of columns in\n",
    "matrix A is equal to the number of rows in matrix B.\n",
    "\n",
    "2. Secondly, the first row of matrix A multiplied by the first column of matrix B\n",
    "gives us the first element in the matrix AB, and so on.\n",
    "\n",
    "3. Thirdly, when multiplying, order matters — specifically, AB ≠ BA.\n",
    "\n",
    "4. Lastly, the element at row i, column j is the product of the ith row of matrix A and the jth column of matrix B.\n",
    "\n",
    "The order in which we multiply matters. We must keep the matrices\n",
    "in order, but we do have some flexibility. As we can see in the following equation, the parentheses can be moved:\n",
    "\n",
    "![](https://imgur.com/7GL8N5q.png)\n",
    "\n",
    "\n",
    "### Implementing a dot production with numpy\n",
    "\n",
    "#### this will create integer as array-elements\n",
    "\n",
    "```py\n",
    "import numpy as np\n",
    "\n",
    "A = np.random.randint(10, _size_=(2,3))\n",
    "B = np.random.randint(10, _size_=(3,2))\n",
    "\n",
    "This will create 2 Matrices as below e.g.\n",
    "[[9 0 8] [3 5 5]]\n",
    "\n",
    "[[5 2]\n",
    " [0 7]\n",
    " [0 7]\n",
    "]\n",
    "\n",
    "print(np.dot(A, B))\n",
    "\n",
    "# Output\n",
    "[[45 74]\n",
    "[15 76]]\n",
    "\n",
    "```\n",
    "#### The Hadamard Product — Other Kinds of Matrix Multiplication\n",
    "\n",
    "Hadamard multiplication is defined for matrices of the same shape as the multiplication of each element of one matrix by the corresponding element of the other matrix. Hadamard multiplication is often denoted by as below, for two matrices A(n×m) and B(n×m) we have\n",
    "\n",
    "![](https://imgur.com/OCHZD2v.png)\n",
    "\n",
    "\n",
    "**So in general Mathematic form for the single independent variable case**\n",
    "\n",
    "![](https://imgur.com/mPTHGZ4.png)\n",
    "\n",
    "#### **Multiple Linear Regression (MLR)**\n",
    "\n",
    "Suppose that the response variable Y and at least one predictor variable xi are quantitative. Then the equation for a specific Y value under the **MLR** model is\n",
    "\n",
    "![](https://imgur.com/UcQBTKp.png)\n",
    "\n",
    "\n",
    "**And now in matrix notation, these n sets of equations become**\n",
    "\n",
    "![](https://cdn-images-1.medium.com/max/800/1*gytWOZW1m99KNjzv5dCNeg.png)\n",
    "\n",
    "\n",
    "where **Y** is the vector of the response variable and is an n × 1 vector of dependent variables, **X** is the matrix of the k independent/explanatory variables (usually the first column is a column of ones for the constant term) and is an n × p matrix of predictors, **β** is a p × 1 vector of unknown coefficients, and e is an n × 1 vector of unknown errors. Equivalently,\n",
    "\n",
    "![](https://imgur.com/WcoZ69u.png)\n",
    "\n",
    "\n",
    "Where\n",
    "**Y** : is output vector for n **training examples.\n",
    "X** : is matrix of size n***p** where each ith row belongs to ith training set.\n",
    "**β** : is weight vector of size p for p training features.\n",
    "\n",
    "### Different types of Matrix Differentiation\n",
    "\n",
    "#### **1-Differentiation with Respect to a Scalar**\n",
    "\n",
    "Differentiation of a structure (vector or matrix, for example) with respect to a scalar is quite simple; it just yields the ordinary derivative of each element of the structure in the same structure. **Thus, the derivative of a vector or a matrix with respect to a scalar variable is a vector or a matrix, respectively, of the derivatives of the individual elements**.\n",
    "\n",
    "**1.A-Derivatives of Vectors with Respect to Scalars**\n",
    "\n",
    "The derivative of the vector y(x) = (y1 , . . . , yn ) with respect to the scalar x\n",
    "is the vector\n",
    "\n",
    "![](https://imgur.com/yDkanT8.png)\n",
    "\n",
    "\n",
    "The second or higher derivative of a vector with respect to a scalar is likewise a vector of the derivatives of the individual elements; that is, it is an array of higher rank.\n",
    "\n",
    "**1.B-Derivatives of Matrices with Respect to Scalars**\n",
    "\n",
    "The derivative of the matrix **Y(x)** defined as below\n",
    "\n",
    "![](https://imgur.com/uDsUz9Q.png)\n",
    "\n",
    "\n",
    "The second or higher derivative of a matrix with respect to a scalar is\n",
    "likewise a matrix of the derivatives of the individual elements.\n",
    "\n",
    "#### 1.C — Derivatives of Functions with Respect to Scalars\n",
    "\n",
    "Differentiation of a function of a vector or matrix that is linear in the elements\n",
    "of the vector or matrix involves just the differentiation of the elements, fol-\n",
    "lowed by application of the function. For example, the derivative of a trace of\n",
    "a matrix is just the trace of the derivative of the matrix. On the other hand,\n",
    "the derivative of the determinant of a matrix is not the determinant of the\n",
    "derivative of the matrix\n",
    "\n",
    "#### 1.D — Higher-Order Derivatives with Respect to Scalars\n",
    "\n",
    "Because differentiation with respect to a scalar does not change the rank of the object (“rank” here means rank of an array or “shape”), higher-order derivatives\n",
    "\n",
    "![](https://imgur.com/mEEP8Ey.png)\n",
    "\n",
    "\n",
    "if those derivatives exist. This vector is called the gradient of the scalar-valued\n",
    "function, and is sometimes denoted by ∇f (x)\n",
    "\n",
    "#### 2.B — Derivatives of Vectors with Respect to Vectors; The Jacobian\n",
    "\n",
    "The derivative of an m-vector-valued function of an n-vector argument consists of nm scalar derivatives. These derivatives could be put into various structures. Two obvious structures are an n × m matrix and an m × n matrix.\n",
    "\n",
    "For a function,\n",
    "\n",
    "![](https://imgur.com/Z6HcNaA.png)\n",
    "\n",
    "\n",
    "if those derivatives exist. This derivative is called the matrix gradient and\n",
    "is denoted by ∇f for the vector-valued function f . (Note that the ∇\n",
    "symbol can denote either a vector or a matrix, depending on whether the\n",
    "function being differentiated is scalar-valued or vector-valued.)\n",
    "\n",
    "**The m × n matrix is called the Jacobian of f and is denoted by Jf as below**\n",
    "\n",
    "![](https://imgur.com/MZmwDGD.png)\n",
    "\n",
    "\n",
    "So we have m = n functions and parameters, in this case. Generally speaking, though, **the Jacobian matrix is the collection of all m × n possible partial derivatives (m rows and n columns), which is the stack of m gradients with respect to x:**\n",
    "\n",
    "![](https://imgur.com/oacOqf6.png)\n",
    "\n",
    "\n",
    "So if we are predicting house-price with the above MLR equation, then _θ_0 will be the basic/base price of a house, then _θ_1 as the price per room, _θ_2 as the price per KM-distance from the nearest Airport.\n",
    "\n",
    "Using the definition of matrix multiplication, our multivariate hypothesis function can be concisely represented as:\n",
    "\n",
    "![](https://imgur.com/Y5y5PUp.png)\n",
    "\n",
    "\n",
    "So now the Cost function takes the following form\n",
    "\n",
    "![](https://imgur.com/lbPyKhu.png)\n",
    "\n",
    "Where the thetas θ are the weights, and the above partial derivative for any weights wj will be as below\n",
    "\n",
    "![](https://imgur.com/x97JW0Q.png)\n",
    "\n",
    "\n",
    "So the Gradient-Descent process for Multivariate case becomes\n",
    "\n",
    "![](https://imgur.com/C6pxz05.png)\n",
    "\n",
    "\n",
    "And that's why we take the transpose of θ to multiply with column-vector x to get the hypothesis (as earlier mentioned in this article)\n",
    "\n",
    "![](https://imgur.com/MT7a9f9.png)\n",
    "\n",
    "\n",
    "### Refresher — Matrix-Derivative Identities required for the Mathematical Derivation of the Gradient of a Matrix w.r.t. to Vectors\n",
    "\n",
    "The derivative of a matrix is usually referred to as the gradient and denoted as ∇. Consider a function\n",
    "\n",
    "![](https://imgur.com/3j5xQwK.png)\n",
    "\n",
    "![](https://imgur.com/7om526N.png)\n",
    "\n",
    "---\n",
    "\n",
    "### Matrix Transpose Identities\n",
    "\n",
    "![](https://imgur.com/7DlM2cH.png)\n",
    "\n",
    "\n",
    "Because of the associativity of matrix multiplication, this relation can be\n",
    "extended as\n",
    "\n",
    "![](https://imgur.com/3yGUHVn.png)\n",
    "\n",
    "\n",
    "### Now proceed to find the Gradient of the Cost function.\n",
    "\n",
    "#### First a Refresher on the Gradient Descent Algorithm\n",
    "\n",
    "To implement Gradient Descent, you need to compute the gradient of the cost function with regard to each model parameter θj. In other words, you need to calculate how much the cost function will change if you change θj just a little bit. This is called a partial derivative. It is like asking “What is the slope of the mountain under my feet if I face east?”. and then asking the same question facing north.\n",
    "\n",
    "Now you would recognize the very well-known cost function\n",
    "\n",
    "![](https://imgur.com/6suRuNj.png)\n",
    "\n",
    "\n",
    "Notice that this formula involves calculations over the full training set **X**, at each Gradient Descent step! This is why the algorithm is called **Batch Gradient Descent**: it uses the whole batch of training data at every step.\n",
    "\n",
    "Once you have the gradient vector, which points uphill, just go in the opposite direction to go downhill. This means subtracting ∇θMSE(θ) from θ. This is where the learning rate η comes into play:5 multiply the gradient vector by η to determine the size of the downhill step\n",
    "\n",
    "![](https://imgur.com/ASVJX8Y.png)\n",
    "\n",
    "\n",
    "Where\n",
    "**Y** : is output vector for n **training examples.\n",
    "X** : is matrix of size n***p** where each ith row belongs to ith training set.\n",
    "**β** : is weight vector of size p for p training features.\n",
    "\n",
    "Note that β in the above is not a scalar, but a vector.\n",
    "\n",
    "Now we have the RSS defined as\n",
    "\n",
    "![](https://imgur.com/kScgoJE.png)\n",
    "\n",
    "\n",
    "Then for the following assumption\n",
    "\n",
    "![](https://imgur.com/0Zrl67y.png)\n",
    "\n",
    "\n",
    "Then for the whole matrix (i.e. the whole set of training data set or the whole set of Hypothesis Equation ), we will get\n",
    "\n",
    "![](https://imgur.com/KI0GW8J.png)\n",
    "\n",
    "\n",
    "Which in the final Vector form\n",
    "\n",
    "![](https://imgur.com/RIFnFp4.png)\n",
    "\n",
    "Note, that in the last equality, I had to get the Transpose of **X** because when doing matrix multiplication — that's a dot product of rows of the first matrix to columns of the second matrix. The number of columns of the 1st matrix must equal the number of rows of the 2nd matrix.\n",
    "\n",
    "So by transposing the _p-_th column of **X** ends up being the _p-_th row of the X-Transposed. Thus, when doing a dot product between the p-th row of **X-Transposed with (y — Xβ) it will match perfectly as I am** using all of the entries of the _p-_th column of **X**\n",
    "\n",
    "### Alternative-2 — And below is an alternative calculation for arriving at the same Vectorized Gradient formulae for training-data in Matrix form.\n",
    "\n",
    "Here, I am denoting the coefficients with **θ or Theta (instead of β that we used above in our Alternative-1 Gradient Calculation — only to make the presentation differentiable)**\n",
    "\n",
    "Again assume we have our Training Set of data as below\n",
    "\n",
    "![](https://imgur.com/Pv7ha4Z.png)\n",
    "\n",
    "So we can say the below\n",
    "\n",
    "And now again, we need to use the same **vector identity** mentioned above, that for a vector z, we have\n",
    "\n",
    "![](https://imgur.com/RZHXy6H.png)\n",
    "\n",
    "\n",
    "We have already introduced the trace operator of a Matrix, written as **“tr.”** Now we need to use a couple of more **matrix derivatives Identities** (that I am just stating below here, and they all have robust Mathematical proofs, the details of which I am not including here).\n",
    "\n",
    "So below 2 M**atrix Derivative Identities hold true and we need to use them to arrive at the Gradient Calculation.**\n",
    "\n",
    "![](https://imgur.com/cCx8lWo.png)\n",
    "\n",
    "\n",
    "So now Final Gradient Calculation will be as below\n",
    "\n",
    "![](https://imgur.com/OdD10oL.png)\n",
    "\n",
    "\n",
    "And also the below Matrix Identity\n",
    "\n",
    "![](https://imgur.com/Dl3pwMU.png)\n",
    "\n",
    "And above is the exact formulae that we will implement in Python/Numpy very soon below.\n",
    "\n",
    "### **So now let's go back to the original Cost Function**\n",
    "\n",
    "![](https://imgur.com/thWUKY1.png)\n",
    "\n",
    "Compare the above with the Gradient-Descent formulae for the Numerical case\n",
    "\n",
    "![](https://imgur.com/hUylz9H.png)\n",
    "\n",
    "\n",
    "### A manual example of the Gradient-Descent implementation\n",
    "\n",
    "Let's say for simple single variable training dataset we have the following values\n",
    "\n",
    "```\n",
    "x,y1,12,23,34,4\n",
    "```\n",
    "\n",
    "Further, assume,\n",
    "\n",
    "α (learning rate) = 1\n",
    "\n",
    "m (number of training examples) =4\n",
    "\n",
    "Setting _θ_0 to 0 and _θ_1 to 1\n",
    "\n",
    "So we have the linear equation\n",
    "\n",
    "![](https://imgur.com/LEDwj0M.png)\n",
    "\n",
    "\n",
    "at the _i_−th row\n",
    "\n",
    "Further I denote,\n",
    "\n",
    "![](https://imgur.com/EOlJlCw.png)\n",
    "\n",
    "\n",
    "So, regardless of how many times I apply the GD algorithm, the value of θ1 will be constantly equal to 1, since at every iteration we have _θ_0=0 and _θ_1=1\n",
    "\n",
    "#### Now extend the above to a multivariable case,\n",
    "\n",
    "Assume theta values have been picked at random as below\n",
    "\n",
    "![](https://imgur.com/klTA8vc.png)\n",
    "\n",
    "\n",
    "And the GD-Algorithm is,\n",
    "\n",
    "![](https://imgur.com/STnCpdd.png)\n",
    "\n",
    "---\n",
    "\n",
    "### Python/Numpy Implementation\n",
    "\n",
    "Please refer to the jupyter notebook\n",
    "\n",
    "First, **generate a training dataset** in Matrix form\n",
    "\n",
    "\n",
    "**NumPy zeros()** function in above **—** you can create an array that only contains only zeros using the NumPy zeros() function with a specific shape. The shape is row by column format. Its syntax is as below\n",
    "\n",
    "![](https://imgur.com/HK87s4H.png)\n",
    "\n",
    "\n",
    "For example, the code to generate a Matrix of 2 by 3 (2 rows and 3 columns)\n",
    "\n",
    "import numpy as np\n",
    "\n",
    "A = np.zeros(shape = (2, 3))\n",
    "print\n",
    "\n",
    "# Output below\n",
    "[[0. 0. 0.]\n",
    " [0. 0. 0.]]\n",
    "\n",
    "Which produces an array like the following:\n",
    "\n",
    "![](https://imgur.com/YuDc7ca.png)\n",
    "\n",
    "\n",
    "If I run the above **gen_data()** function above for a set of 5 training data-set as below with bias and variance of 20 and 10 respectively\n",
    "\n",
    "gen_data(5, 20, 10)\n",
    "\n",
    "I will have the following form of output\n",
    "\n",
    "```\n",
    "(array([[1., 0.],\n",
    " [1., 1.],\n",
    " [1., 2.],\n",
    " [1., 3.],\n",
    " [1., 4.]]),\n",
    " array([22.38023816, 24.18406356, 28.01360908, 26.80051617, 29.30101971])\n",
    " )\n",
    "\n",
    "```\n",
    "And now the function for Gradient-Descent implementing the Grdient formulae for a Mactrix that we derived above\n",
    "\n",
    "![](https://imgur.com/p3U0Ab6.png)\n",
    "\n",
    "\n",
    "**The** [**full notebook is here**](https://gist.github.com/rohan-paul/de3999b34b480d5ee26acc9105bf50a9)\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "#### Reference and Further Reading\n",
    "\n",
    "**Matrix Multiplication** — [https://en.wikipedia.org/wiki/Matrix_multiplication](https://en.wikipedia.org/wiki/Matrix_multiplication)\n",
    "\n",
    "**Matrix-Calculus —** [https://en.wikipedia.org/wiki/Matrix_calculus](https://en.wikipedia.org/wiki/Matrix_calculus)\n",
    "\n",
    "**Vector_Field** — [https://en.wikipedia.org/wiki/Vector_field](https://en.wikipedia.org/wiki/Vector_field)\n",
    "\n",
    "**Matrix Transpose Properties** -[https://en.wikipedia.org/wiki/Transpose#Properties](https://en.wikipedia.org/wiki/Transpose#Properties)\n",
    "\n",
    "**Matrix Cookbook —** [https://www.math.uwaterloo.ca/~hwolkowi/matrixcookbook.pdf](https://www.math.uwaterloo.ca/~hwolkowi/matrixcookbook.pdf)\n",
    "\n",
    "**Online Calculation of Matrix Derivative** — [http://www.matrixcalculus.org/](http://www.matrixcalculus.org/)"
   ]
  }
 ],
 "metadata": {
  "language_info": {
   "name": "python"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}