This notebook presents how to manipulate matrices in NumPy. NumPy is the fundamental package for scientific computing with Python.
import numpy as np
NumPy provides the type N-dimensional array.
a = np.array([[1, 2], [3, 4]])
type(a)
numpy.ndarray
a
array([[1, 2], [3, 4]])
You can multiply two numpy.ndarray
s:
a * a
array([[ 1, 4], [ 9, 16]])
But you can see that arrays are not exactly mathematical matrices. The expected result for matrix multiplication would be np.array([[7, 10], [15, 22]])
. Operators act element-wise on arrays. Since this is the case in matrix addition too, you can use arrays as matrices when you do:
a + a
array([[2, 4], [6, 8]])
Adding an array and a scalar follows the expected (sensible) behaviour.
a + 1
array([[2, 3], [4, 5]])
NumPy provides the usual functions to apply on arrays.
np.sum(a) # Sums all elements
10
Partial sum along the first (index 0) axis:
np.sum(a, 0) # | (axis=0)
array([4, 6])
Partial sum along the second (index 1) axis:
np.sum(a, 1) # --> (axis=1)
array([3, 7])
np.sum(a, axis=1)
array([3, 7])
Exercise
What does np.sum(a, 2)
do?
Likewise, we can call the mean
function:
np.mean(a)
2.5
Over the entire array or partially:
np.mean(a, axis=0)
array([ 2., 3.])
You can look for values in your array:
np.where(a == 2)
(array([0]), array([1]))
returns the coordinates of element 2. Note that coordinates are consistent with the mathematical matrix convention. Coordinates in this format can be readily used to access elements; trivially:
coord = np.where(a == 1)
a[coord]
array([1])
Note that everything is stored in arrays. If you have several instances, it looks like the following:
a1 = np.ones((2, 2))
print "Array is", a1
coord1 = np.where(a1 == 1)
print "Ones are found at (row, column) =", coord1
print "Values at these coordinates", a1[coord1]
print "Values of former array at these coordinates", a[coord1]
Array is [[ 1. 1.] [ 1. 1.]] Ones are found at (row, column) = (array([0, 0, 1, 1]), array([0, 1, 0, 1])) Values at these coordinates [ 1. 1. 1. 1.] Values of former array at these coordinates [1 2 3 4]
If you have no instances, say:
np.where(a == 0)
(array([], dtype=int64), array([], dtype=int64))
There is no element equal to 0 in a
, so we get an empty array.
Exercise
Access the element of a
which lies on the second row and first column.
To select the entire second row, use a slice as follows:
a[1, :]
array([3, 4])
Slices illustrate the power of high-level programming. Do not write for loops. Let the computer worry about how to do the element-by-element operations!
Matrix multiplication can be expressed by the dot product of 2D arrays:
np.dot(a, a)
array([[ 7, 10], [15, 22]])
The dot product of 1D arrays expresses inner product of vectors:
b = np.array([1, 2, 3])
np.dot(b, b) # 1*1 + 2*2 + 3*3
14
Of course, dot()
only works on arrays with compatible shapes. For example,
np.dot(a, b)
--------------------------------------------------------------------------- ValueError Traceback (most recent call last) <ipython-input-22-579c274cec9b> in <module>() ----> 1 np.dot(a, b) ValueError: objects are not aligned
a.shape
(2, 2)
b.shape
(3,)
You cannot multiply a 2-by-2 matrix with a vector of size 3.
Exercise
Create an object named c
such that you can compute np.dot(a, c)
.
Try np.dot(c, a)
. What do you notice?
# NumPy does not distinguish between row and colum vectors.
If you want to use *
as the operator for matrix multiplication, you need to create matrix objects.
m = np.matrix([[1, 2], [3, 4]])
m * m
matrix([[ 7, 10], [15, 22]])
You can convert an array into a matrix and vice versa.
am = np.matrix(a)
ma = np.array(m)
Matrix has all the features of array. You want to use the matrix type if your problem is linear algebra. Indeed, vectors are then 1-by-N matrices.
bm = np.matrix(b)
bm.shape
(1, 3)
Transpose it to get a column vector:
np.transpose(bm)
matrix([[1], [2], [3]])
Otherwise, typically if you are representing multi-dimensional grids, you should use array.
Always remember to look at the docs: http://docs.scipy.org/doc/numpy/
You will find implementations for conjugate
, convolve
, correlate
, diagonal
, fft
, gradient
, ... These functions are faster than anything you could easily write and ... someone else has tested and debugged them! :)