Matrices#

Import#

import numpy as np

Matrix as a NumPy array#

Everything in NumPy is an array. A matrix is also an array. Let us create a simple matrix:

\[\begin{split} \textbf{M} = \begin{bmatrix} 1 & 2 & 3\\ 4 & 5 & 6\\ 7 & 8 & 9 \end{bmatrix} \end{split}\]

In NumPy:

M = np.array([
    [1, 2, 3],
    [4, 5, 6],
    [7, 8, 9]
])
M
array([[1, 2, 3],
       [4, 5, 6],
       [7, 8, 9]])

Adding two matrices#

Let us now add the following matrices:

\[\begin{split} \textbf{A} = \begin{bmatrix} 1 & 2\\ 3 & 4 \end{bmatrix}, \textbf{B} = \begin{bmatrix} 5 & 6\\ 7 & 8 \end{bmatrix} \end{split}\]

then,

\[\begin{split} \textbf{C} = \textbf{A} + \textbf{B} = \begin{bmatrix} 6 & 8\\ 10 & 12 \end{bmatrix} \end{split}\]

In NumPy:

A = np.array([
    [1, 2],
    [3, 4]])
B = np.array([
    [5, 6],
    [7, 8]
])
C = A + B
C
array([[ 6,  8],
       [10, 12]])

Scaling a matrix#

Scaling a matrix is nothing but element-wise multiplication:

\[\begin{split} \textbf{M} = \begin{bmatrix} 1 & 2 & 3\\ 4 & 5 & 6\\ 7 & 8 & 9 \end{bmatrix} \end{split}\]

then,

\[\begin{split} 3 \textbf{M} = \begin{bmatrix} 3 & 6 & 9\\ 12 & 15 & 18\\ 21 & 24 & 27 \end{bmatrix} \end{split}\]

In NumPy:

M = np.array([
    [1, 2, 3],
    [4, 5, 6],
    [7, 8, 9]
])
3 * M
array([[ 3,  6,  9],
       [12, 15, 18],
       [21, 24, 27]])

Element-wise multiplication of matrices#

Consider two matrices:

\[\begin{split} \textbf{A} = \begin{bmatrix} 1 & 2\\ 3 & 4 \end{bmatrix}, \textbf{B} = \begin{bmatrix} 5 & 6\\ 7 & 8 \end{bmatrix} \end{split}\]

The element-wise product is given by \(\textbf{A} \odot \textbf{B}\):

\[\begin{split} \textbf{C} = \textbf{A} \odot \textbf{B} = \begin{bmatrix} 5 & 12\\ 21 & 32 \end{bmatrix} \end{split}\]

In NumPy:

A = np.array([
    [1, 2],
    [3, 4]
])
B = np.array([
    [5, 6],
    [7, 8]
])
C = A * B
C
array([[ 5, 12],
       [21, 32]])

Element-wise functions of matrices#

Given a matrix, we sometimes would want to apply a function to every element of the matrix. We will consider two examples.

Example-1#

For example, we may want to take the absolute value of all the elements. Let us say \(f(x) = |x|\), then:

\[\begin{split} \mathbf{A} = \begin{bmatrix} -1 & 2\\ -3 & -4 \end{bmatrix} \end{split}\]

then:

\[\begin{split} \begin{bmatrix} f(-1) & f(2)\\ f(-3) & f(-4) \end{bmatrix} = \begin{bmatrix} 1 & 2\\ 3 & 4 \end{bmatrix} \end{split}\]

In NumPy, this becomes:

A = np.array([
    [-1, 2],
    [-3, -4]
])
np.abs(A)
array([[1, 2],
       [3, 4]])

Example-2#

We might want to square each element of the matrix. If \(\textbf{A}\) is a matrix, then \(\textbf{B}\) could be defined element-wise as follows:

\[ B_{ij} = A_{ij}^2 \]

Let us compute \(\mathbf{B}\) for the following matrix:

\[\begin{split} \mathbf{A} = \begin{bmatrix} 1 & \sqrt{2}\\ \sqrt{3} & 2 \end{bmatrix} \end{split}\]

In NumPy:

A = np.array([
    [1, np.sqrt(2)],
    [np.sqrt(3), 2]
])
A ** 2
array([[1., 2.],
       [3., 4.]])

Transpose of a matrix#

Given a matrix \(\textbf{M}\):

\[\begin{split} \textbf{M} = \begin{bmatrix} 1 & 2 & 3\\ 4 & 5 & 6 \end{bmatrix} \end{split}\]

then, its transpose \(\textbf{M}^{T}\) is:

\[\begin{split} \textbf{M}^{T} = \begin{bmatrix} 1 & 4\\ 2 & 5\\ 3 & 6 \end{bmatrix} \end{split}\]

In NumPy:

M = np.array([
    [1, 2, 3],
    [4, 5, 6]
])
np.transpose(M)
array([[1, 4],
       [2, 5],
       [3, 6]])

Alternatively:

M = np.array([
    [1, 2, 3],
    [4, 5, 6]
])
M.transpose()
array([[1, 4],
       [2, 5],
       [3, 6]])

Equivalently, we have:

M.T
array([[1, 4],
       [2, 5],
       [3, 6]])

This is what we will be sticking to given that it is very close to the corresponding math notation.

Shape and dimension of a matrix#

Matrices are “two dimensional” arrays. So all matrices in NumPy have array-dimension equal to two. The shape of the NumPy array gives what we usually call the dimension of the matrix in the linear algebra sense.

Explore these two ideas for:

\[\begin{split} \mathbf{M} = \begin{bmatrix} 1 & 2 & 3\\ 4 & 5 & 6\\ \end{bmatrix} \end{split}\]

In NumPy:

M = np.array([
    [1, 2, 3],
    [4, 5, 6]
])
print(M.shape)
print(M.ndim)
(2, 3)
2

Vectors as matrices#

Each vector can be viewed as a matrix. Column vectors are matrices of shape \((d, 1)\). Row vectors are matrices of shape \((1, d)\). Let us look at how NumPy treats both these cases:

x = np.array([
    [1],
    [2],
     [3]
])
print(x.shape)
print(x.ndim)
x
(3, 1)
2
array([[1],
       [2],
       [3]])
x = np.array([
    [1, 2, 3]
])
print(x.shape)
print(x.ndim)
x
(1, 3)
2
array([[1, 2, 3]])

To create a row or column vector from a 1D NumPy array, we can use np.newaxis. Focus on the syntax for now; the reason we are doing this will become clear when we discuss indexing and slicing.

# Column vector
x = np.array([1, 2, 3])
x = x[:, np.newaxis]
x.shape
(3, 1)
# Row vector
x = np.array([1, 2, 3])
x = x[np.newaxis, :]
x.shape
(1, 3)

Alternatively, we can also use a method called reshape:

# Column vector
x = np.array([1, 2, 3])
x = x.reshape(3, 1)
x.shape
(3, 1)
# Row vector
x = np.array([1, 2, 3])
x = x.reshape(1, 3)
x.shape
(1, 3)

Products involving matrices and vectors#

We will look at the following products:

  • matrix - matrix

  • matrix - vector

  • vector - matrix

  • vector - vector

Product of two matrices#

Given two matrices:

\[\begin{split} \textbf{A} = \begin{bmatrix} 1 & 2 & 3\\ 4 & 5 & 6 \end{bmatrix}, \textbf{B} = \begin{bmatrix} 6 & 7\\ 8 & 9\\ 10 & 11 \end{bmatrix} \end{split}\]

then,

\[\begin{split} \textbf{C} = \textbf{A} \times \textbf{B} = \begin{bmatrix} 52 & 58\\ 124 & 139 \end{bmatrix} \end{split}\]

In NumPy:

A = np.array([
    [1, 2, 3],
    [4, 5, 6]
])
B = np.array([
    [6, 7],
    [8, 9],
    [10, 11]
])
C = A @ B
C
array([[ 52,  58],
       [124, 139]])

Product of a matrix and a (column) vector#

Given the matrix \(\mathbf{A}\) and the vector \(\mathbf{x}\):

\[\begin{split} \mathbf{A} = \begin{bmatrix} 1 & 2 & 3\\ 4 & 5 & 6\\ 7 & 8 & 9 \end{bmatrix}, \mathbf{x} = \begin{bmatrix} 6\\ 7\\ 8 \end{bmatrix} \end{split}\]

The product \(\mathbf{Ax}\) is given by:

\[\begin{split} \mathbf{C} = \mathbf{A x} = \begin{bmatrix} 44\\ 107\\ 170 \end{bmatrix} \end{split}\]

In NumPy:

A = np.array([
    [1, 2, 3],
    [4, 5, 6],
    [7, 8, 9]
])
x = np.array([6, 7, 8])
C = A @ x
C
array([ 44, 107, 170])

Product of a (row) vector and a matrix#

Given the matrix \(\mathbf{A}\) and the vector \(\mathbf{x}\):

\[\begin{split} \mathbf{A} = \begin{bmatrix} 1 & 2 & 3\\ 4 & 5 & 6\\ 7 & 8 & 9 \end{bmatrix}, \mathbf{x} = \begin{bmatrix} 6\\ 7\\ 8 \end{bmatrix} \end{split}\]

The product \(\mathbf{x}^T \mathbf{A}\) is given by:

\[ \mathbf{C} = \mathbf{x}^T \mathbf{A} = \begin{bmatrix} 90 & 111 & 132 \end{bmatrix} \]

In NumPy:

A = np.array([
    [1, 2, 3],
    [4, 5, 6],
    [7, 8, 9]
])
x = np.array([6, 7, 8])
x @ A
array([ 90, 111, 132])

(Inner) Product of a (row) vector and a (column) vector#

The product of a row vector and a column vector is nothing but the usual dot product:

\[\begin{split} \mathbf{x}^T = \begin{bmatrix} 1 & 2 & 3 \end{bmatrix}, \quad \mathbf{y} = \begin{bmatrix} 4\\ 5\\ 6 \end{bmatrix} \end{split}\]

The product \(\mathbf{x}^T \mathbf{y}\) is then:

\[ \mathbf{x}^T \mathbf{y} = 32 \]

In NumPy:

x = np.array([1, 2, 3])
y = np.array([4, 5, 6])
x @ y
32

(Outer) Product of a (column) vector and a (row) vector#

The product of a column vector and a row vector is an outer product:

\[\begin{split} \mathbf{x} = \begin{bmatrix} 1\\ 2\\ 3 \end{bmatrix}, \quad \mathbf{y} = \begin{bmatrix} 4 & 5 & 6 \end{bmatrix} \end{split}\]

The product \(\mathbf{x} \mathbf{y}^T\) is then:

\[\begin{split} \mathbf{x} \mathbf{y}^T = \begin{bmatrix} 4 & 5 & 6\\ 8 & 10 & 12\\ 12 & 15 & 18 \end{bmatrix} \end{split}\]

In NumPy:

x = np.array([1, 2, 3])
y = np.array([4, 5, 6])
np.outer(x, y)
array([[ 4,  5,  6],
       [ 8, 10, 12],
       [12, 15, 18]])

Equivalently:

x = np.array([1, 2, 3])
y = np.array([4, 5, 6])
x = x[:, np.newaxis]
y = y[np.newaxis, :]
x @ y
array([[ 4,  5,  6],
       [ 8, 10, 12],
       [12, 15, 18]])

Matrix of zeros#

In many algorithms, we might have to initialize a matrix with zeros. For example, consider a \(2 \times 4\) matrix:

\[\begin{split} \mathbf{M} = \begin{bmatrix} 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 \end{bmatrix} \end{split}\]

In NumPy:

np.zeros((2, 4))
array([[0., 0., 0., 0.],
       [0., 0., 0., 0.]])

Matrix of ones#

Similar to a matrix of zeros, we can come up with a matrix of ones.

\[\begin{split} \mathbf{M} = \begin{bmatrix} 1 & 1\\ 1 & 1\\ 1 & 1 \end{bmatrix} \end{split}\]

In NumPy:

np.ones((3, 2))
array([[1., 1.],
       [1., 1.],
       [1., 1.]])

Identity matrix#

Often, we might have to deal with identity matrices. A \(3 \times 3\) identity matrix is as follows:

\[\begin{split} \mathbf{I} = \begin{bmatrix} 1 & 0 & 0\\ 0 & 1 & 0\\ 0 & 0 & 1 \end{bmatrix} \end{split}\]

In NumPy:

np.eye(5)
array([[1., 0., 0., 0., 0.],
       [0., 1., 0., 0., 0.],
       [0., 0., 1., 0., 0.],
       [0., 0., 0., 1., 0.],
       [0., 0., 0., 0., 1.]])

Diagonal matrices#

Another special kind of matrix. Let us create the following matrix:

\[\begin{split} \mathbf{D} = \begin{bmatrix} 1 & 0 & 0 & 0\\ 0 & 2 & 0 & 0\\ 0 & 0 & 3 & 0\\ 0 & 0 & 0 & 4 \end{bmatrix} \end{split}\]

In NumPy:

np.diag([1, 2, 3, 4])
array([[1, 0, 0, 0],
       [0, 2, 0, 0],
       [0, 0, 3, 0],
       [0, 0, 0, 4]])

Indexing and Slicing#

Just like lists in Python, NumPy arrays can be indexed and sliced. Slicing is useful if we want to work with a portion of an array. We will look at some examples.

Example-1: Indexing#

Consider:

\[\begin{split} M = \begin{bmatrix} 1 & 3 & 0\\ 4 & 2 & 5\\ 9 & 8 & 7 \end{bmatrix} \end{split}\]

To get the element \(5\), we can do the following:

M = np.array([
    [1, 3, 0],
    [4, 2, 5],
    [9, 8, 7]
])
M[1][2]
5

Alternatively, we can also use:

M = np.array([
    [1, 3, 0],
    [4, 2, 5],
    [9, 8, 7]
])
M[1, 2]
5

Note that NumPy uses zero-indexing, just like Python. In general, M[row, col] is the element in the row row and columnn col assuming zero indexing.

Example-2: Row-slice#

We will extract the third row of the matrix \(\mathbf{M}\):

\[\begin{split} \mathbf{M} = \begin{bmatrix} 1 & 2\\ 3 & 4\\ 5 & 6\\ 7 & 8\\ 9 & 10 \end{bmatrix} \end{split}\]

In NumPy:

M = np.arange(1, 11).reshape(5, 2)
M
array([[ 1,  2],
       [ 3,  4],
       [ 5,  6],
       [ 7,  8],
       [ 9, 10]])
M[2, :]
array([5, 6])

The syntax is to be understood as follows. To extract the third row, we fix the row index to be \(2\) (zero-indexing) and get hold of all elements in that row. To do this, we use a : in the place of the column index so that it allows all elements to come in.

Example-3: Column slice#

Let us now extract the second column of the following matrix:

\[\begin{split} \mathbf{M} = \begin{bmatrix} 1 & 2 & 3\\ 4 & 5 & 6\\ 7 & 8 & 9 \end{bmatrix} \end{split}\]

In NumPy:

M = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
M[:, 1]
array([2, 5, 8])

Here, we want the second column. So we fix the column index to \(1\) and set the row index to :, indicating that all elements in that column should be included.

Example-4: Submatrix slice#

Now, we want to extract the \(2 \times 2\) submatrix colored in blue from \(M\):

\[\begin{split} \mathbf{M} = \begin{bmatrix} 1 & 2 & 3 & 4\\ 5 & 6 & \color{blue}7 & \color{blue}8\\ 9 & 10 & \color{blue}{11} & \color{blue}{12}\\ 13 & 14 & 15 & 16 \end{bmatrix} \end{split}\]

In NumPy:

M = np.array([
    [1, 2, 3, 4],
    [5, 6, 7, 8],
    [9, 10, 11, 12],
    [13, 14, 15, 16]
])
M[1:3, 2:4]
array([[ 7,  8],
       [11, 12]])

Start point of the slice is included and the end-point is excluded, just as we do in native Python.

Other Matrix Operations#

There are several important matrix operations that we will list down here. The np.linalg module helps us perform these operations effortlessly.

Rank#

We can get the rank of a matrix as follows. For example:

\[\begin{split} \mathbf{M} = \begin{bmatrix} 1 & 2 & 3 & 4\\ 5 & 6 & 7 & 8\\ 9 & 10 & 11 & 12\\ 13 & 14 & 15 & 16 \end{bmatrix} \end{split}\]
M = np.array([
    [1, 2, 3, 4],
    [5, 6, 7, 8],
    [9, 10, 11, 12],
    [13, 14, 15, 16]
])
np.linalg.matrix_rank(M)
2

Inverse#

We can get the inverse as follows. For example:

\[\begin{split} M = \begin{bmatrix} 3 & -4\\ 4 & 3 \end{bmatrix} \end{split}\]
M = np.array([[3, -4], [4, 3]])
np.linalg.inv(M)
array([[ 0.12,  0.16],
       [-0.16,  0.12]])

Pseuodoinverse#

The pseudoinverse \(\mathbf{M}^{\dagger}\) of a matrix \(\mathbf{M}\) can be computed using NumPy.

\[\begin{split} \mathbf{M} = \begin{bmatrix} 1 & 2 & 3\\ 3 & 6 & 9 \end{bmatrix} \end{split}\]
M = np.array([[1, 2, 3], [3, 6, 9]])
np.linalg.pinv(M)
array([[0.00714286, 0.02142857],
       [0.01428571, 0.04285714],
       [0.02142857, 0.06428571]])

The pseudoinverse reduces to the inverse for square invertible matrices.

Eigenvalues and Eigenvectors#

Given a symmetric matrix \(\mathbf{M}\) let us find its eigenvalues and eigenvectors:

\[\begin{split} \mathbf{M} = \begin{bmatrix} 1 & 0 & -3\\ 0 & 5 & 2\\ -3 & 2 & 8 \end{bmatrix} \end{split}\]
M = np.array([[1, 0, -3], [0, 5, 2], [-3, 2, 8]])
eigval, eigvec = np.linalg.eigh(M)
eigval
array([-0.20942046,  4.36588492,  9.84353554])

The eigenvalues are arranged in ascending order. The corresponding eigenvectors appear in the columns of the matriix eigvec. The eigenpairs are given below:

eigval[0], eigvec[:, 0]
eigval[1], eigvec[:, 1]
eigval[2], eigvec[:, 2];