Numerical programming with numpy

Read then Launch

This content is best viewed in html because jupyter notebook cannot display some content (e.g. figures, equations) properly. You should finish reading this page first and then launch it as an interactive notebook in Google Colab (faster, Google account needed) or Binder by clicking the rocket symbol () at the top.

In python, we can use the numpy package to perform basic matrix operations. The numpy package is a fundamental package for scientific computing with Python. It contains among other things:

  • A powerful \(N\)-dimensional array object.

  • Sophisticated (broadcasting) functions.

  • A collection of routines for linear algebra, Fourier transform, and random number generation.

  • A collection of routines for numerical integration and optimisation.

  • Key tools for working with numerical data in Python.

  • Support for large, multi-dimensional arrays and matrices.

Creating matrices/arrays

Create a matrix using list

# create matrices
x = [[1, 2, 3], [4, 5, 6]]
print(x)
[[1, 2, 3], [4, 5, 6]]

Create an array using numpy.array

import numpy as np

# create an 1D array

x = np.array([1, 2, 3, 4, 5])
print(x)
[1 2 3 4 5]

The shape attribute of an array object returns a tuple consisting of array dimensions. For a matrix with n rows and m columns, its shape will be (n,m). The length of the shape tuple is therefore the number of axes, ndim.

print(x.shape)
(5,)

np.arange() is used to create a sequence of numbers. It is similar to the built-in range() function, but returns an ndarray instead of a list. The arguments are (start, stop, step), which are the same as for range(), but it also accepts float arguments. If step is not given, it defaults to 1. If start is not given, it defaults to 0. If stop is not given, it defaults to start and start is set to 0. Note that np.arange() excludes the right end of range specification.

x = np.arange(1, 11)
print(x)
[ 1  2  3  4  5  6  7  8  9 10]
# note: np.arange actually can result in unexpected results; check np.arange(0.2, 0.6, 0.4) vs np.arange(0.2, 1.6, 1.4)
print(np.arange(0.2, 0.6, 0.4))
print(np.arange(0.2, 1.6, 1.4))
[0.2]
[0.2 1.6]

The reshape() method can be used to reshape an array.

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

x = np.reshape(x, (2, 3))
print(x)
[1 2 3 4 5 6]
[[1 2 3]
 [4 5 6]]

Try to run the following cell multiple times, can you get the same output every time?

# generate random numbers/array/matrices

mu, sigma = 0, 1
a = np.random.randint(0, 10, 5)
b = np.random.random((2, 3))

print(a)
print(b)

x = np.random.normal(mu, sigma, 5)
y = x + np.random.normal(20, 0.1, 5)

print(x)
print(y)
[9 2 0 4 6]
[[0.30045579 0.34532917 0.60851618]
 [0.73175115 0.02666414 0.01931966]]
[-1.4779433  -1.33043671 -1.54065478 -1.61139813 -0.98688384]
[18.51214766 18.57095824 18.57120978 18.4464537  19.09870361]

Run the following cell multiple times, can you get the same output every time?

np.random.seed(123)

print(np.random.normal(mu, sigma, 5))
[-1.0856306   0.99734545  0.2829785  -1.50629471 -0.57860025]

Basic matrix operations

  • Transpose a matrix (array).

The T attribute or transpose() method can be used to transpose an array.

a = np.array([[1, 2], [3, 4]])

print(a)

print(a.T)

print(a.transpose())
[[1 2]
 [3 4]]
[[1 3]
 [2 4]]
[[1 3]
 [2 4]]
  • Add a scalar to an array, or multiply an array by a scalar.

print(a + 1)
print(a - 1)
print(a * 2)
print(a / 2)
[[2 3]
 [4 5]]
[[0 1]
 [2 3]]
[[2 4]
 [6 8]]
[[0.5 1. ]
 [1.5 2. ]]
  • Add or multiply two arrays of the same size element-wise.

b = np.array([[5, 6], [7, 8]])

print(a + b)

print(a * b)

print(np.multiply(a, b))
[[ 6  8]
 [10 12]]
[[ 5 12]
 [21 32]]
[[ 5 12]
 [21 32]]
  • Matrix multiplication.

# product of two matrices
c = np.array([[1, 2, 3], [4, 5, 6]])

print(np.dot(a, c))
[[ 9 12 15]
 [19 26 33]]
#  product of matrix and vector
d = np.array([1, 2])

print(np.dot(a, d))
[ 5 11]

Statistics with NumPy

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

print(np.sqrt(x))

print(x**2)

print(np.square(x))
[1.         1.41421356 1.73205081 2.         2.23606798 2.44948974]
[ 1  4  9 16 25 36]
[ 1  4  9 16 25 36]
x = np.random.randint(0, 10, 5)

print("The mean value of the array x is: ", np.mean(x))
print("The median value of the array x is: ", np.median(x))
print("The standard deviation of the array x is: ", np.std(x))
print("The variance of the array x is: ", np.var(x))
print("The max value of the array x is: ", np.min(x))
print("The min value of the array x is: ", np.max(x))
The mean value of the array x is:  3.4
The median value of the array x is:  1.0
The standard deviation of the array x is:  3.49857113690718
The variance of the array x is:  12.239999999999998
The max value of the array x is:  0
The min value of the array x is:  9

Indexing in NumPy

A = np.arange(1, 17, 1).reshape(4, 4).transpose()
print(A)
[[ 1  5  9 13]
 [ 2  6 10 14]
 [ 3  7 11 15]
 [ 4  8 12 16]]

Get an element from a matrix, for example, get the fourth element in the third row (index starts from 0).

print(A[2, 3])
15
# to select a submatrix, need the non-singleton dimension of your indexing array to be aligned with the axis you're indexing into,
# e.g. for an n x m 2D subarray: A[n by 1 array,1 by m array]
A[[[0], [2]], [1, 3]]
array([[ 5, 13],
       [ 7, 15]])

Using the : operator can get a row or a column.

# the last two examples include either no index for the columns or no index for the rows. These indicate that Python should include all columns or all rows, respectively
A[0, :]
array([ 1,  5,  9, 13])

The default indexing is row-wise, and therefore the : for column can be omitted.

A[0]
array([ 1,  5,  9, 13])

Using start:stop;step_size to get a submatrix.

  • start and stop are the indices of the first and last elements of the submatrix, and step_size is the step size.

  • start and stop can be negative, which means the indices are counted from the right.

  • step_size can be negative, which means the submatrix is reversed.

  • start and stop can be omitted, which means the submatrix starts from the first element or ends at the last element.

  • step_size can be omitted, which means the step size is 1.

# this is another way to do that
A[0:3:2, 1:4:2]
array([[ 5, 13],
       [ 7, 15]])
# select all columns in those two rows
A[0:3:2, :]
array([[ 1,  5,  9, 13],
       [ 3,  7, 11, 15]])
# select all row in those two columns
A[:, 1:4:2]
array([[ 5, 13],
       [ 6, 14],
       [ 7, 15],
       [ 8, 16]])
# '-' sign has a different meaning and good usage in Python. This means index from the end, -1 means the last element
A[-1, -1]
16

Boolean indexing. There are other ways to let Python keep all rows except certain indices.

# index using a boolean value
ind = np.ones((4,), bool)
ind[[0, 2]] = False
print(ind)
[False  True False  True]
print(A[ind, :])

print(A[ind])
[[ 2  6 10 14]
 [ 4  8 12 16]]
[[ 2  6 10 14]
 [ 4  8 12 16]]

Exercises

1. \(A\) = \(\begin{bmatrix} 8 & 5 & 9 \\ 7 & 19 & 14 \end{bmatrix}\); \(B\) = \( \begin{bmatrix} 6 & 13 & 15 \\ 17 & 22 & 12 \\ 0 & 3 & 9 \end{bmatrix}\)

i. \(AB\) = ?

ii. \(2A\) = ?

iii. \(AB + 2A\) = ?

You have to do it programmatically. Match the result with Linear Algebra Excercise 2.

# Write your code below to answer the question

Compare your answer with the reference solution below

import numpy as np

a = np.array([[8, 5, 9], [7, 19, 14]])
b = np.array([[6, 13, 15], [17, 22, 12], [0, 3, 9]])

print("i.")
print(np.dot(a, b))

print("ii.")
print(2 * a)

print("iii.")
print(np.dot(a, b) + (2 * a))
i.
[[133 241 261]
 [365 551 459]]
ii.
[[16 10 18]
 [14 38 28]]
iii.
[[149 251 279]
 [379 589 487]]

2. Create two \(3 \times 3\) matrices using the NumPy random() function. Output the result of multiplying the \(2\) matrices (as a NumPy array). Don’t forget to use the seed() function to create reproducible results.

# Write your code below to answer the question

Compare your answer with the reference solution below

import numpy as np

np.random.seed(123)
x = np.random.random((3, 3))
y = np.random.random((3, 3))

print(np.dot(x, y))
[[0.56600587 0.29748762 0.66145826]
 [0.84396992 0.30934893 0.76255113]
 [1.03984767 0.46521201 1.07199874]]