import numpy as np
import time
First, a motivating example for why we might want a NUMerical PYthon (numpy) library.
num_elements = 1000000
input_array = np.arange(num_elements)
start_time = time.time()
return_array = [0] * len(input_array)
for key, val in enumerate(input_array):
return_array[key] = val * val
print(time.time() - start_time)
0.21277904510498047
start_time = time.time()
return_array_vectorized = np.power(input_array, 2)
print(time.time() - start_time)
0.0031232833862304688
Array operations
Array construction
np.zeros(5)
array([0., 0., 0., 0., 0.])
np.identity(3)
array([[1., 0., 0.], [0., 1., 0.], [0., 0., 1.]])
array_1 = np.array([[1,2],[3,4],[5,6]]) # An array with arbitrary elements.
array_1
array([[1, 2], [3, 4], [5, 6]])
Array access
array_1[1,1] # Access an element
4
array_1[:,1] # Access a row
array([2, 4, 6])
Array broadcasting of addition and scalar multiplication
array_2 = np.add(array_1, 1) # array_1 + 1
array_2
array([[2, 3], [4, 5], [6, 7]])
np.multiply(array_1, 2) # array_1 * 2
array([[ 2, 4], [ 6, 8], [10, 12]])
array_1 / 2.0 # np.multiply(array_1, 1/2.0)
array([[0.5, 1. ], [1.5, 2. ], [2.5, 3. ]])
Array elementwise addition
array_1 + array_2 # np.add(array_1, array_2)
array([[ 3, 5], [ 7, 9], [11, 13]])
Array elementwise multiplication
np.multiply(array_1, array_1) # array_1 * array_1
array([[ 1, 4], [ 9, 16], [25, 36]])
Dot-product
array_3 = np.array([[1, 2]])
array_4 = np.array([[3], [4]])
np.dot(array_3, array_4) # Inner-product
array([[11]])
array_5 = np.dot(array_4, array_3) # Outer-product
array_5
array([[3, 6], [4, 8]])
Array info
array_5.sum(axis=0)
array([ 7, 14])
array_5.mean(axis=0)
array([3.5, 7. ])
array_5.std(axis=0)
array([0.5, 1. ])
array_5.max(axis=0)
array([4, 8])
array_5.shape
(2, 2)
Transpose
M1 = array_1.T
M1
array([[1, 3, 5], [2, 4, 6]])
Matrix elementwise exponentiation
np.power(M1, 2) # M1 ^ 2?
array([[ 1, 9, 25], [ 4, 16, 36]])
Matrix multiplication
np.matmul(M1.T, M1) # or np.dot(M1.T, M1)
array([[ 5, 11, 17], [11, 25, 39], [17, 39, 61]])
Matrix multiplication is not commutative
print(np.matmul(M1, M1.T))
print(np.matmul(M1.T, M1))
[[35 44] [44 56]] [[ 5 11 17] [11 25 39] [17 39 61]]
Matrix Inverse
np.linalg.inv(np.array([[2, 0], [0, 2]]))
array([[0.5, 0. ], [0. , 0.5]])
Inverse fails if matrix is singular.
np.linalg.inv(np.matmul(M1.T, M1)) # This should raise an exception.
--------------------------------------------------------------------------- LinAlgError Traceback (most recent call last) <ipython-input-65-9e5bbb853a28> in <module>() ----> 1 np.linalg.inv(np.matmul(M1.T, M1)) # This should raise an exception. ~/miniconda3/envs/py37/lib/python3.7/site-packages/numpy/linalg/linalg.py in inv(a) 530 signature = 'D->D' if isComplexType(t) else 'd->d' 531 extobj = get_linalg_error_extobj(_raise_linalgerror_singular) --> 532 ainv = _umath_linalg.inv(a, signature=signature, extobj=extobj) 533 return wrap(ainv.astype(result_t, copy=False)) 534 ~/miniconda3/envs/py37/lib/python3.7/site-packages/numpy/linalg/linalg.py in _raise_linalgerror_singular(err, flag) 87 88 def _raise_linalgerror_singular(err, flag): ---> 89 raise LinAlgError("Singular matrix") 90 91 def _raise_linalgerror_nonposdef(err, flag): LinAlgError: Singular matrix
Evaluating linear systems Ax = b for b
A1 = np.array([[1, 1], [0, 1]])
x1 = np.array([[2], [2]])
b1 = np.matmul(A1, x1)
b1
array([[4], [2]])
Solving linear systems Ax = b for x
x1_verify = np.linalg.solve(A1, b1)
all(x1 == x1_verify)
True
If A is singular then the linear system may be overdetermined (no solution) when solving for x.
A2 = np.array([[1, 1], [2, 2]]) # Note that the rank of this matrix is 1.
x2 = np.array([[2], [3]])
potential_sol = np.matmul(A2, x2)
potential_sol
array([[ 5], [10]])
b2 = np.array([[5], [11]])
x2_verify = np.linalg.solve(A2, b2) # This should raise an exception.
# There is no selection of x to solve this linear system for A2 and b2.
--------------------------------------------------------------------------- TypeError Traceback (most recent call last) <ipython-input-69-afbd29c7fcd1> in <module>() 1 b2 = np.array([[5], 11]) ----> 2 x2_verify = np.linalg.solve(A2, b2) # This should raise an exception. 3 # There is no selection of x to solve this linear system for A2 and b2. ~/miniconda3/envs/py37/lib/python3.7/site-packages/numpy/linalg/linalg.py in solve(a, b) 392 signature = 'DD->D' if isComplexType(t) else 'dd->d' 393 extobj = get_linalg_error_extobj(_raise_linalgerror_singular) --> 394 r = gufunc(a, b, signature=signature, extobj=extobj) 395 396 return wrap(r.astype(result_t, copy=False)) TypeError: No loop matching the specified signature and casting was found for ufunc solve1
x2_verify = np.linalg.solve(A2, potential_sol) # This should raise an exception.
# There happens to exist an x to solve this system, but numpy still fails.
--------------------------------------------------------------------------- LinAlgError Traceback (most recent call last) <ipython-input-70-ee72a78afdea> in <module>() ----> 1 x2_verify = np.linalg.solve(A2, potential_sol) # This should raise an exception. 2 # There happens to exist an x to solve this system, but numpy still fails. ~/miniconda3/envs/py37/lib/python3.7/site-packages/numpy/linalg/linalg.py in solve(a, b) 392 signature = 'DD->D' if isComplexType(t) else 'dd->d' 393 extobj = get_linalg_error_extobj(_raise_linalgerror_singular) --> 394 r = gufunc(a, b, signature=signature, extobj=extobj) 395 396 return wrap(r.astype(result_t, copy=False)) ~/miniconda3/envs/py37/lib/python3.7/site-packages/numpy/linalg/linalg.py in _raise_linalgerror_singular(err, flag) 87 88 def _raise_linalgerror_singular(err, flag): ---> 89 raise LinAlgError("Singular matrix") 90 91 def _raise_linalgerror_nonposdef(err, flag): LinAlgError: Singular matrix
Another linear system that is underdetermined (many solutions) when solving for x
A3 = np.array([[1, 0, 0], [0, 1, 1]])
x3 = np.array([[1], [1], [1]])
b3 = np.matmul(A3, x3)
b3
array([[1], [2]])
np.linalg.solve(A3, b3) # This should raise an exception.
--------------------------------------------------------------------------- LinAlgError Traceback (most recent call last) <ipython-input-72-e0ec3f51dfe0> in <module>() ----> 1 np.linalg.solve(A3, b3) # This should raise an exception. ~/miniconda3/envs/py37/lib/python3.7/site-packages/numpy/linalg/linalg.py in solve(a, b) 379 a, _ = _makearray(a) 380 _assertRankAtLeast2(a) --> 381 _assertNdSquareness(a) 382 b, wrap = _makearray(b) 383 t, result_t = _commonType(a, b) ~/miniconda3/envs/py37/lib/python3.7/site-packages/numpy/linalg/linalg.py in _assertNdSquareness(*arrays) 213 m, n = a.shape[-2:] 214 if m != n: --> 215 raise LinAlgError('Last 2 dimensions of the array must be square') 216 217 def _assertFinite(*arrays): LinAlgError: Last 2 dimensions of the array must be square