Tutorial_06b_MatrixAlgebra.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Basic Matrix Algebra\n",
    "\n",
    "This notebook presents some basic linear algebra in Julia."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "printyellow (generic function with 1 method)"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "using Printf, LinearAlgebra\n",
    "\n",
    "include(\"jlFiles/printmat.jl\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Adding and Multiplying: A Matrix and a Scalar"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With a matrix $A$ and a scalar $c$, do\n",
    "\n",
    "1. `A*c` (textbook: $Ac$) to multiply each element of $A$ by $c$\n",
    "\n",
    "2. `A .+ c` (textbook: $A+cJ$, where $J$ is a matrix of ones) to add $c$ to each element of $A$, and similarly `A .- c` ($A-c$)\n",
    "\n",
    "Watch out when the number comes first: `2.+A` is not allowed since it is ambiguous. However, `2.0.+A` and `2 .+ A` both work."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "A:\n",
      "     1         3    \n",
      "     3         4    \n",
      "\n",
      "c:\n",
      "    10    \n",
      "\n",
      "A*c:\n",
      "    10        30    \n",
      "    30        40    \n",
      "\n",
      "A .+ c:\n",
      "    11        13    \n",
      "    13        14    \n",
      "\n"
     ]
    }
   ],
   "source": [
    "A = [1 3;3 4]\n",
    "c = 10\n",
    "\n",
    "println(\"A:\")\n",
    "printmat(A)\n",
    "println(\"c:\")\n",
    "printmat(c)\n",
    "\n",
    "println(\"A*c:\")\n",
    "printmat(A*c)\n",
    "\n",
    "println(\"A .+ c:\")\n",
    "printmat(A .+ c)          #notice the dot in .+"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Adding and Multiplying Two Matrices"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With two matrices of the same dimensions ($A$ and $B$), do\n",
    "\n",
    "`A+B` (textbook: $A+B$) to add them (element by element), and similarly `A-B` (textbook: $A-B$).\n",
    "\n",
    "Multiplying matrices ($A$ and $B$) of conformable dimensions\n",
    "\n",
    "`A*B` (textbook: $AB$)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "A:\n",
      "     1         3    \n",
      "     3         4    \n",
      "\n",
      "B:\n",
      "     1         2    \n",
      "     3        -2    \n",
      "\n",
      "A+B:\n",
      "     2         5    \n",
      "     6         2    \n",
      "\n",
      "A*B:\n",
      "    10        -4    \n",
      "    15        -2    \n",
      "\n"
     ]
    }
   ],
   "source": [
    "A = [1 3;3 4]               #A and B are 2x2 matrices\n",
    "B = [1 2;3 -2]\n",
    "println(\"A:\")\n",
    "printmat(A)\n",
    "println(\"B:\")\n",
    "printmat(B)\n",
    "\n",
    "println(\"A+B:\")\n",
    "printmat(A+B)\n",
    "\n",
    "println(\"A*B:\")\n",
    "printmat(A*B)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Transpose\n",
    "\n",
    "You can transpose a numerical matrix `A` by `A'`. \n",
    "\n",
    "Notice that (in Julia) `A` and `B=A'` share the same elements (changing one changes the other). If you want an independent copy, use `B=copy(A')`.\n",
    "\n",
    "For an array of other elements (for instance, strings), use `permutedims(A)` to swap the dimensions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "A: \n",
      "     1         2         3    \n",
      "     4         5         6    \n",
      "\n",
      "A': \n",
      "     1         4    \n",
      "     2         5    \n",
      "     3         6    \n",
      "\n"
     ]
    }
   ],
   "source": [
    "A = [1 2 3;4 5 6]\n",
    "println(\"A: \")\n",
    "printmat(A)\n",
    "println(\"A': \")\n",
    "printmat(A')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Matrix Inverse"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A matrix inverse of an $nxn$ matrix $A$:\n",
    "\n",
    "`inv(A)` or `A^(-1)` (textbook: $A^{-1}$)\n",
    "\n",
    "The inverse is such that $AA^{-1}=I$ and $A^{-1}A=I$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "A:\n",
      "     1         3    \n",
      "     3         4    \n",
      "\n",
      "inv(A):\n",
      "    -0.800     0.600\n",
      "     0.600    -0.200\n",
      "\n",
      "inv(A)*A:\n",
      "     1.000    -0.000\n",
      "     0.000     1.000\n",
      "\n"
     ]
    }
   ],
   "source": [
    "A = [1 3;3 4]\n",
    "println(\"A:\")\n",
    "printmat(A)\n",
    "\n",
    "println(\"inv(A):\")\n",
    "printmat(inv(A))\n",
    "\n",
    "println(\"inv(A)*A:\")\n",
    "printmat(inv(A)*A)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## The Identity Matrix\n",
    "\n",
    "The identity matrix $I_n$ can often be represented by `I` and then Julia will compare with the surrounding code to create the right dimension. For instance, if `A` is a square matrix, then `I+A` works.\n",
    "\n",
    "If you still need to specify the dimension, then `1I(3)` or `Matrix(1I(3))` will create an $I_3$ matrix. (The former is an explicitly diagonal matrix, will the second is a full matrix.) In Julia 1.6 in will be possible to do `1I[1:3,1:3]`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "I + A\n",
      "     2         3    \n",
      "     3         5    \n",
      "\n",
      "Matrix(1I(3))\n",
      "     1         0         0    \n",
      "     0         1         0    \n",
      "     0         0         1    \n",
      "\n"
     ]
    }
   ],
   "source": [
    "println(\"I + A\")\n",
    "printmat(I + A)\n",
    "\n",
    "println(\"Matrix(1I(3))\")\n",
    "printmat(Matrix(1I(3)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Vectors: Inner and Outer Products\n",
    "\n",
    "There are several different ways to think about a vector in mathematics: as a $K \\times 1$ matrix (a column vector), a $1 \\times K$ matrix (a row vector) or just a flat $K-$vector. Julia uses flat vectors but they are mostly interchangable with column vectors. \n",
    "\n",
    "The inner product of two (column) vectors with $k$ elements is calculated as `x'z` or `dot(x,y)` (textbook: $x'z$ or $x \\cdot z$) to get a scalar. (You can also use  or `x⋅y` where the dot is obtained by `\\cdot + TAB`, but that is sometimes hard to distinguish from  or `x.y`) \n",
    "\n",
    "In contrast, the outer of two (column) vectors with $k$ elements is calculated\n",
    "(textbook: $xz'$) to get a $k\\times k$ matrix."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x and z\n",
      "    10         2    \n",
      "    11         5    \n",
      "\n",
      "x'z: \n",
      "    75    \n",
      "x*z':\n",
      "    20        50    \n",
      "    22        55    \n",
      "\n"
     ]
    }
   ],
   "source": [
    "x = [10,11]                  #[10;11] gives the same\n",
    "z = [2,5]\n",
    "println(\"x and z\")\n",
    "printmat([x z])\n",
    "\n",
    "println(\"x'z: \")\n",
    "printlnPs(x'z)               #dot(x,z) gives the same\n",
    "\n",
    "println(\"x*z':\")\n",
    "printmat(x*z')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Vectors: Quadratic Forms\n",
    "\n",
    "A quadratic form ($A$ is an $n \\times n$ matrix and x is an $n$ vector): `x'A*x` (textbook: $x'Ax$) to get a scalar. \n",
    "\n",
    "(In Julia 1.5+ there is also the form `dot(x,A,x)`)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x:\n",
      "    10    \n",
      "    11    \n",
      "\n",
      "A:\n",
      "     1         3    \n",
      "     3         4    \n",
      "\n",
      "x'A*x: \n",
      "  1244    \n"
     ]
    }
   ],
   "source": [
    "A = [1 3;3 4]\n",
    "x = [10,11]\n",
    "\n",
    "println(\"x:\")\n",
    "printmat(x)\n",
    "println(\"A:\")\n",
    "printmat(A)\n",
    "\n",
    "println(\"x'A*x: \")\n",
    "printlnPs(x'A*x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Vectors: Extracting Vectors from Matrices\n",
    "\n",
    "Notice that `A[1,:]` and `A[:,1]` both give flat vectors. In case you want a row vector use `A[1:1,:]`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2-element Vector{Int64}:\n",
       " 1\n",
       " 3"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "A[1,:]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## OLS Notation\n",
    "\n",
    "$X'X$ or $\\sum\\nolimits_{t=1}^{T}x_{t}x_{t}^{\\prime}$?\n",
    "\n",
    "Let $x_t$ be a (column) vector with values of $K$ regressors for observation $t$. Then $x_{t}x_{t}^{\\prime}$ is the outer product (a $K\\times K$ matrix) and $\\sum\\nolimits_{t=1}^{T}x_{t}x_{t}^{\\prime}$ is just the sum (of each element) across the $T$ observations.\n",
    "\n",
    "We can calculate the same thing by (a) letting $X$ be a $T\\times K$ matrix with $x_{t}^{\\prime}$ in row $t$ and (b) then do $X'X$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "X\n",
      "     1.000    -1.000\n",
      "     1.000     0.000\n",
      "     1.000     1.000\n",
      "\n",
      "sum of outer products, three versions\n",
      "     3.000     0.000\n",
      "     0.000     2.000\n",
      "\n",
      "     3.000     0.000\n",
      "     0.000     2.000\n",
      "\n",
      "     3.000     0.000\n",
      "     0.000     2.000\n",
      "\n"
     ]
    }
   ],
   "source": [
    "x₁ = [1,-1]     #a (column) vector\n",
    "x₂ = [1,0]\n",
    "x₃ = [1,1.0]\n",
    "\n",
    "X = [x₁';x₂';x₃']\n",
    "\n",
    "println(\"X\")\n",
    "printmat(X)\n",
    "\n",
    "(T,K) = (size(X,1),size(X,2))\n",
    "\n",
    "Sxx1 = x₁*x₁' + x₂*x₂' + x₃*x₃'      #just to illustrate\n",
    "\n",
    "Sxx2 = zeros(K,K)\n",
    "for t = 1:T\n",
    "    #global Sxx2                   #only needed in script\n",
    "    Sxx2 = Sxx2 + X[t,:]*X[t,:]'   #X[t,:] becomes a flat vector\n",
    "end\n",
    "\n",
    "Sxx3 = X'X\n",
    "\n",
    "println(\"sum of outer products, three versions\")\n",
    "printmat(Sxx1)\n",
    "printmat(Sxx2)\n",
    "printmat(Sxx3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "@webio": {
   "lastCommId": null,
   "lastKernelId": null
  },
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Julia 1.6.0",
   "language": "julia",
   "name": "julia-1.6"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "1.6.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}