In this chapter we work with vectors and matrices. We cover the basic calculations that can be done with both. We also discuss how vectors and matrices can be manipulated.
A vector is a list of numbers that we can do “math” with. The numbers in the vector are indexed, so that we can access them. Note that vector indexing in Python starts with zero, not with 1. So Python counts the first element in a vector as element 0, the second as element 1, etc. We first need to import some important packages so that Python understands basic numerical procedures and definitions, like vectors and matrices and not with a list. The command np.array from the numpy package does just that.
Here are some examples. Vectors x y year and names are assigned as follows:
import numpy as np
# You can document your script files using the ``#`` symbol.
# This allows you to add commentary to your codes.
xv = np.array([1, 3.3333, 4.8976, 9])
yv = np.array([9.233, 0.3453, 3.29999, 2.8])
print("xv = {}".format(xv))
print("yv = {}".format(yv))
print("xv[0] = {}".format(xv[0]))
print("xv[1] = {}".format(xv[1]))
print("xv[2] = {}".format(xv[2]))xv = [1. 3.3333 4.8976 9. ]
yv = [9.233 0.3453 3.29999 2.8 ]
xv[0] = 1.0
xv[1] = 3.3333
xv[2] = 4.8976# You can document your script files using the ``#`` symbol.
# This allows you to add commentary to your codes.
xv <- c(1, 3.3333, 4.8976, 9)
yv <- c(9.233, 0.3453, 3.29999, 2.8)
cat("xv =", xv, "\n")
cat("yv =", yv, "\n")
cat("xv[1] =", xv[1], "\n")
cat("xv[2] =", xv[2], "\n")
cat("xv[3] =", xv[3], "\n")xv = 1 3.3333 4.8976 9
yv = 9.233 0.3453 3.29999 2.8
xv[1] = 1
xv[2] = 3.3333
xv[3] = 4.8976 Remember that Python starts numbering the elements with 0 and not 1! So if you want to access the first element in a vector x you’d need to type: xv[0].
If you would like to print the vector formatted, you need to set the format option of the numpy object as:
Note that if you like a vector: 2000,2001,2002,2003 you need to write the arange command with the upper limit+1.
If you would like to generate a certain number of elements between two boundary points you can use np.linspace(). Let us say we would like to generate some numbers between the lower bound of 2000 and the upper bound of 2003. We call this spanning a grid from 2000 to 2003 where we basically break the interval into evenly spaced sub-intervals. This is important for graphing (or plotting) mathematical functions which we do in the next chapter.
Let us say we want to span a grid with 5 gridpoints, so that the interval from 2000 to 2003 is broken in equally wide sub-intervals, we would write the following:
Remember from the previous chapter how arrays are very similar to lists and tuples. An array is basically a list of numbers, and numbers only! This means that indexing will work on an array, just as it worked on the list object. Same with slicing if you remember that from the previous chapter.
Here is again a list example. Let us define it first and then extract some objects from the list (which is itself an object – remember in Python everything is an object).
# A list of strings is simply:
names_list = ["Tom", "Dick", "Harry", "Patrick"]
# A tuple of strings is simply:
names_tuple = ("Tom", "Dick", "Harry", "Patrick")
print("names_list= {}".format(names_list))
print("names_tuple= {}".format(names_tuple))
print("names_list[1]= {}".format(names_list[1]))
print("names_list[3]= {}".format(names_list[3]))names_list= ['Tom', 'Dick', 'Harry', 'Patrick']
names_tuple= ('Tom', 'Dick', 'Harry', 'Patrick')
names_list[1]= Dick
names_list[3]= Patrick# A list of strings is simply:
names_list <- c("Tom", "Dick", "Harry", "Patrick")
# A vector of strings is equivalent to a tuple in Python
names_vector <- c("Tom", "Dick", "Harry", "Patrick")
cat("names_list =", names_list, "\n")
cat("names_vector =", names_vector, "\n")
cat("names_list[2] =", names_list[2], "\n")
cat("names_list[4] =", names_list[4], "\n")names_list = Tom Dick Harry Patrick
names_vector = Tom Dick Harry Patrick
names_list[2] = Dick
names_list[4] = Patrick Let us do some simple math with arrays (or vectors, which are one dimensional arrays).
import numpy as np
# Element-by-element operations
x1 = np.array([1,3,4,9])
x2 = np.array([2,5,6,3])
#
print(" --- OUTPUT: --- ")
print("x1= {}".format(x1))
print("x2= {}".format(x2))
print(" ---------------- ")
print("x1+x2= {}".format(x1+x2))
print("x1*x2= {}".format(x1*x2)) --- OUTPUT: ---
x1= [1 3 4 9]
x2= [2 5 6 3]
----------------
x1+x2= [ 3 8 10 12]
x1*x2= [ 2 15 24 27]# Element-by-element operations
x1 <- c(1, 3, 4, 9)
x2 <- c(2, 5, 6, 3)
cat(" --- OUTPUT: --- \n")
cat("x1 =", x1, "\n")
cat("x2 =", x2, "\n")
cat(" ---------------- \n")
cat("x1 + x2 =", x1 + x2, "\n")
cat("x1 * x2 =", x1 * x2, "\n") --- OUTPUT: ---
x1 = 1 3 4 9
x2 = 2 5 6 3
----------------
x1 + x2 = 3 8 10 12
x1 * x2 = 2 15 24 27 As you can see, the numbers are added element-by-element and they are multiplied element-by-element.
This is different behavior from what we have seen so far from lists. Let us repeat the above example, but instead of numpy-arrays we use lists instead. We would then have:
x1_list = [1,3,4,9]
x2_list = [2,5,6,3]
#
print(" --- OUTPUT: --- ")
print("x1_list = {}".format(x1_list ))
print("x2_list = {}".format(x2_list ))
print(" ---------------- ")
print("x1_list + x2_list = {}".format(x1_list +x2_list )) --- OUTPUT: ---
x1_list = [1, 3, 4, 9]
x2_list = [2, 5, 6, 3]
----------------
x1_list + x2_list = [1, 3, 4, 9, 2, 5, 6, 3]x1_list <- c(1, 3, 4, 9)
x2_list <- c(2, 5, 6, 3)
cat(" --- OUTPUT: --- \n")
cat("x1_list =", x1_list, "\n")
cat("x2_list =", x2_list, "\n")
cat(" ---------------- \n")
cat("x1_list + x2_list =", x1_list + x2_list, "\n") --- OUTPUT: ---
x1_list = 1 3 4 9
x2_list = 2 5 6 3
----------------
x1_list + x2_list = 3 8 10 12 This just joins the lists together. It does not add the numbers in the list.
If, on the other hand, we try to multiply the lists we get:
print("x1_list*x2_list= {}".format(x1_list*x2_list))Error: TypeError: can't multiply sequence by non-int of type 'list'cat("x1_list * x2_list =", x1_list * x2_list, "\n")x1_list * x2_list = 2 15 24 27 This results in an error message. You cannot really do “math” on lists because the content of lists can be numbers, words, other weird objects and it is not guaranteed that doing math with these is possible. Hence, Python will not let you do “math” with lists. You can only do “list things” such as adding stuff to it, printing it, extracting info from it, etc.
Matrices are “two dimensional” vectors. Think of a speadsheet with rows and columns that can only be filled with numbers. In Python we define as matrix again with the array() function. Here is a simple Matrix with 2 rows and 2 columns. We call this a 2-by-2 or \(2\times2\) matrix.
import numpy as np
A = np.array([[2,3],[4,5]])
B = np.array([[2,6],[1,3]])
#
print(" --- OUTPUT: --- ")
print("A= \n {}".format(A))
print("B= \n {}".format(B))
print(" ---------------- ")
print("A*B= \n {}".format(A*B))
print("A-B= \n {}".format(A-B))
print(" ---------------- ") --- OUTPUT: ---
A=
[[2 3]
[4 5]]
B=
[[2 6]
[1 3]]
----------------
A*B=
[[ 4 18]
[ 4 15]]
A-B=
[[ 0 -3]
[ 3 2]]
---------------- A <- matrix(c(2, 4, 3, 5), nrow = 2, byrow = TRUE)
B <- matrix(c(2, 1, 6, 3), nrow = 2, byrow = TRUE)
cat(" --- OUTPUT: --- \n")
cat("A= \n", A, "\n")
cat("B= \n", B, "\n")
cat(" ---------------- \n")
cat("A * B= \n", A * B, "\n")
cat("A - B= \n", A - B, "\n")
cat(" ---------------- \n") --- OUTPUT: ---
A=
2 3 4 5
B=
2 6 1 3
----------------
A * B=
4 18 4 15
A - B=
0 -3 3 2
---------------- You can again see that when you add or multiply matrices, the operation is done element-by-element. This means that the number in position 0,0 (row 1, column 1) in matrix A is multiplied with the number in the same position 0,0 (row 1, column 1) in matrix B.
Python has some built in functions to help you generate often used matrices automatically, so you do not always have to type them in.
If you want to generate matrices filled with either zeros or ones of a particular size in Python you can use the following:
import numpy as np
A = np.zeros((3,5), float)
B = np.ones((4,3), float)
print(" --- OUTPUT: --- ")
print("A= \n {}".format(A))
print(" ---------------- ")
print("B= \n {}".format(B)) --- OUTPUT: ---
A=
[[0. 0. 0. 0. 0.]
[0. 0. 0. 0. 0.]
[0. 0. 0. 0. 0.]]
----------------
B=
[[1. 1. 1.]
[1. 1. 1.]
[1. 1. 1.]
[1. 1. 1.]]If you want a matrix with integer numbers, then simply define the matrix as datatype int so that:
import numpy as np
A = np.zeros((3,5), int)
b = np.ones((4,3), int)
print(" --- OUTPUT: --- ")
print("A= \n {}".format(A))
print(" ---------------- ")
print("B= \n {}".format(B)) --- OUTPUT: ---
A=
[[0 0 0 0 0]
[0 0 0 0 0]
[0 0 0 0 0]]
----------------
B=
[[1. 1. 1.]
[1. 1. 1.]
[1. 1. 1.]
[1. 1. 1.]]If you want to create an identity matrix with the value one in the main diagonal and zeros everywhere else you can:
A vector of order \(n>0\) is a set of ordered numbers. \(a = \left[ \begin{array}{c} 4 \\ 3 \end{array} \right], \: e_{1} = \left[ \begin{array}{c} 1 \\ 0 \\ 0 \end{array} \right], \: ...\)
These are column vectors. Row vectors are transposed column vectors, that is \(a'=[4 \: 3], \: e'_{1}=[1 \: 0 \: 0], \: ...\)
In Python vectors are formed as column vectors by default. We can simply transpose vectors using the transpose function t() as follows:
import numpy as np
import matplotlib.pyplot as plt
import math as m
from scipy import stats as st
import time # Imports system time module to time your script
plt.close('all') # close all open figureslibrary(ggplot2) # Load the ggplot2 library for plotting
graphics.off() # Close all open figures
# Other libraries are automatically loaded in R when needed
# There's no direct equivalent to the 'time' module in Python for timing in R.
# However, you can use the 'system.time()' function to measure the execution time of expressions.
# You can proceed with the rest of your R code here.#| results: hold
a = np.array([4,3]) # column vector
aprime = a.transpose() # row vector
print("a= {}".format(a))a= [4 3]print("a'= {}".format(aprime))a'= [4 3]A vector (at least a two dimensional one) has a convenient geometric representation. It is an arrow, where the two coordinates indicate the direction and length of this arrow. Vector \(a = \left[ \begin{array}{c} 4 \\ 3 \end{array} \right]\) points to the upper right (i.e. 4 over, 3 up). The length of a vector can be calculated using the Pythagoras theorem for the triangle. The length of vector \(a = \left[ \begin{array}{c} a_1 \\ a_2 \end{array} \right]\) where \(a_1\) and \(a_2\) are simply numbers, can be calculated as \(\|a\| = \sqrt{a_1^2 + a_2^2}.\) For our example, the vector norm for \(a = \left[ \begin{array}{c} 4 \\ 3 \end{array} \right] \text{ is } \rightarrow \|a\| = \sqrt{4^2 + 3^2} = 5.\) In Python we can simply define a vector and calculate its norm (or length) as
a = np.array([4,3])
norma = np.sqrt(np.sum(a**2))
print("norm(a) = ", norma)
print("norm(a) = ", np.linalg.norm(a)) # built in norm() commandnorm(a) = 5.0
norm(a) = 5.0Adding vectors is simple. We just add all numbers of two vectors ‘element-by-element’. So that \(a + b = \left[ \begin{array}{c} a_1 + b_1 \\ a_2+b_2 \end{array} \right].\) In Python this is done with
a = np.array([4,3])
b = np.array([12,5])
print("a= ", a)
print("b= ", b)
print("a+b=", a+b) # adding up vectors element-by-elementa= [4 3]
b= [12 5]
a+b= [16 8]A vector can be multiplied by a number (we call it a scalar denoted as :math:` lambda` to distinguish it from vectors and the numbers that it contains). The scalar is multiplied with all numbers of the vector. If we multiply the vector with \(\lambda > 1\) then the arrow that the vector symbolizes is becoming longer. If we multiply the vector with \(0 < \lambda < 1\) the arrow gets shorter. If we multiply the vector with a negative number \(\lambda < 0\) then the arrow changes direction. More formally this is \(\lambda * a = \left[ \begin{array}{c} \lambda *a_1 \\ \lambda *a_2 \end{array} \right].\)
In Python it’s simply
a = np.array([4,3])
anorm = np.sqrt(np.sum(a**2)) # norm of vector a
print("a=", a)
print("anorm=", anorm)
b = 5*a
bnorm = np.sqrt(np.sum(b**2)) # norm of longer vector
print("b=", b)
print("bnorm=", bnorm)
c = 0.5*a
cnorm = np.sqrt(np.sum(c**2)) # norm of shorter vector
print("c=", c)
print("cnorm=", cnorm )
d = -3*a # vector changes direction
dnorm = np.sqrt(np.sum(d**2))
print("d=", d)
print("dnorm=", dnorm)a= [4 3]
anorm= 5.0
b= [20 15]
bnorm= 25.0
c= [2. 1.5]
cnorm= 2.5
d= [-12 -9]
dnorm= 15.0a <- c(4, 3)
anorm <- sqrt(sum(a^2)) # norm of vector a
cat("a =", a, "\n")
cat("anorm =", anorm, "\n")
b <- 5 * a
bnorm <- sqrt(sum(b^2)) # norm of longer vector
cat("b =", b, "\n")
cat("bnorm =", bnorm, "\n")
c <- 0.5 * a
cnorm <- sqrt(sum(c^2)) # norm of shorter vector
cat("c =", c, "\n")
cat("cnorm =", cnorm, "\n")
d <- -3 * a # vector changes direction
dnorm <- sqrt(sum(d^2))
cat("d =", d, "\n")
cat("dnorm =", dnorm, "\n")a = 4 3
anorm = 5
b = 20 15
bnorm = 25
c = 2 1.5
cnorm = 2.5
d = -12 -9
dnorm = 15 When multiplying two vectors we form a so called ‘inner product’ as follows: \(a \times b = \sum_{i=1}^{n}a_i*b_i.\)
The following example shows this with numbers:
\[a = \left[ \begin{array}{c} 4 \\ 3 \end{array} \right] \text{ and } b = \left[ \begin{array}{c} 12 \\ 5 \end{array} \right]\]
then
\(a \times b = 4*12 + 3*5.\) In Python the inner product is simply:
A = np.array([[11,12,13,14],[21,22,23,24],[31,32,33,34]])
print("A=", A)
print("---------------")
print()
# Element row 2, column 3 (remember Python starts indexing at 0!!)
print("A[1,2] = ", A[1,2])
print("---------------")
print("A[0,:] = ", A[0,:]) # First row
print("---------------")
print("A[:,0] = ", A[:,0]) # First column
print("---------------")
print("A[1:,:] = ", A[1:,]) # All, except first rowA= [[11 12 13 14]
[21 22 23 24]
[31 32 33 34]]
---------------
A[1,2] = 23
---------------
A[0,:] = [11 12 13 14]
---------------
A[:,0] = [11 21 31]
---------------
A[1:,:] = [[21 22 23 24]
[31 32 33 34]]A <- matrix(c(11, 21, 31, 12, 22, 32, 13, 23, 33, 14, 24, 34), nrow = 3, byrow = TRUE)
cat("A =", A, "\n")
cat("---------------\n\n")
# Element row 2, column 3 (remember R starts indexing at 1!!)
cat("A[2, 3] =", A[2, 3], "\n")
cat("---------------\n")
cat("A[1,] =", A[1,], "\n") # First row
cat("---------------\n")
cat("A[,1] =", A[,1], "\n") # First column
cat("---------------\n")
cat("A[-1,] =", A[-1,], "\n") # All, except first rowA = 11 22 33 21 32 14 31 13 24 12 23 34
---------------
A[2, 3] = 13
---------------
A[1,] = 11 21 31 12
---------------
A[,1] = 11 22 33
---------------
A[-1,] = 22 33 32 14 13 24 23 34 If you want to remove columns or rows you can use the take method and assign the result to a new matrix. Let’s say you would like to eliminate the second column from a matrix A.
print("A=", A)
print("---------------")
# Remove column 2 (or take column 1, 3, and 4)
B = A.take([0,2,3], axis=1)
print("B =", B )A= [[11 12 13 14]
[21 22 23 24]
[31 32 33 34]]
---------------
B = [[11 13 14]
[21 23 24]
[31 33 34]]If you want to replace matrix elements based on certain criteria you can use the criteria as a logical statement in combination with the choose method. In the next example we replace all the elements of the matrix that are larger than 90 with the number 90. We start by setting some elements of the matrix equal to 99.
A[:,0] = 99 # Fill first column with 99
print("A = ", A)
print("---------------")
# Replace all elements > 90 with number 90
B = (A>90).choose(A,90)
print("(A>90).choose(A,90) = ", B)A = [[99 12 13 14]
[99 22 23 24]
[99 32 33 34]]
---------------
(A>90).choose(A,90) = [[90 12 13 14]
[90 22 23 24]
[90 32 33 34]]A[, 1] <- 99 # Fill first column with 99
cat("A =", A, "\n")
cat("---------------\n")
# Replace all elements > 90 with number 90
B <- ifelse(A > 90, 90, A)
cat("(A > 90) * (A, 90) =", B, "\n")A = 99 99 99 21 32 14 31 13 24 12 23 34
---------------
(A > 90) * (A, 90) = 90 90 90 21 32 14 31 13 24 12 23 34 If you want to run a loop “through” a matrix you can first find out the dimensions of the matrix using the shape method and then write nested loops, one for each dimension. Note the indentation of the various print statements in the next example!
A = np.array([[11,12,13,14],[21,22,23,24],[31,32,33,34]])
print("A=", A)
print("---------------")
print(" Start loops ")
print("---------------")
nrRows, nrCols = A.shape
for i in range(nrRows):
print("Row {}".format(i+1))
print("---------------")
for j in range(nrCols):
print("The element in row {} and column {} is: {}".format(i+1, j+1, A[i,j]))
print("Done with row {}".format(i+1))
print("---------------")
print("---------------")
print(" All Done ")
print("---------------")A= [[11 12 13 14]
[21 22 23 24]
[31 32 33 34]]
---------------
Start loops
---------------
Row 1
---------------
The element in row 1 and column 1 is: 11
The element in row 1 and column 2 is: 12
The element in row 1 and column 3 is: 13
The element in row 1 and column 4 is: 14
Done with row 1
---------------
Row 2
---------------
The element in row 2 and column 1 is: 21
The element in row 2 and column 2 is: 22
The element in row 2 and column 3 is: 23
The element in row 2 and column 4 is: 24
Done with row 2
---------------
Row 3
---------------
The element in row 3 and column 1 is: 31
The element in row 3 and column 2 is: 32
The element in row 3 and column 3 is: 33
The element in row 3 and column 4 is: 34
Done with row 3
---------------
---------------
All Done
---------------A <- matrix(c(11, 21, 31, 12, 22, 32, 13, 23, 33, 14, 24, 34), nrow = 3, byrow = TRUE)
cat("A =", A, "\n")
cat("---------------\n")
cat(" Start loops \n")
cat("---------------\n")
nrRows <- nrow(A)
nrCols <- ncol(A)
for (i in 1:nrRows) {
cat("Row ", i, "\n")
cat("---------------\n")
for (j in 1:nrCols) {
cat("The element in row ", i, " and column ", j, " is: ", A[i, j], "\n")
}
cat("Done with row ", i, "\n")
cat("---------------\n")
}
cat("---------------\n")
cat(" All Done \n")
cat("---------------\n")A = 11 22 33 21 32 14 31 13 24 12 23 34
---------------
Start loops
---------------
Row 1
---------------
The element in row 1 and column 1 is: 11
The element in row 1 and column 2 is: 21
The element in row 1 and column 3 is: 31
The element in row 1 and column 4 is: 12
Done with row 1
---------------
Row 2
---------------
The element in row 2 and column 1 is: 22
The element in row 2 and column 2 is: 32
The element in row 2 and column 3 is: 13
The element in row 2 and column 4 is: 23
Done with row 2
---------------
Row 3
---------------
The element in row 3 and column 1 is: 33
The element in row 3 and column 2 is: 14
The element in row 3 and column 3 is: 24
The element in row 3 and column 4 is: 34
Done with row 3
---------------
---------------
All Done
---------------Matrices are “two dimensional” vectors. In Python we define a matrix as
Transposing matrices requires again the t() command. It writes each column as row of a new matrix. So that the transpose of matrix \(A\) in the above example becomes:
If we transpose the transpose of matrix A, we get the original matrix A back.
When adding two matrices \(A\) and \(B\) we simply add all the elements of each matrix ‘element-by-element’. Note that the dimensions (i.e. the number of rows and columns) of the two matrices have to be identical. So if matrix \(A\) has dimension \(m \times n\), that is \(m\) rows and \(n\) columns, then matrix \(B\) needs to be of dimensions \(m \times n\) as well.
A = np.array([[2,3],[4,5]])
B = np.array([[2,6],[1,3]])
C = A+B
print("A=", A)
print("B=", B)
print("C=", C)A= [[2 3]
[4 5]]
B= [[2 6]
[1 3]]
C= [[4 9]
[5 8]]When multiplying two matrices \(A\) and \(B\) we need to make sure that the number of columns of matrix \(A\) is equal the number of rows of matrix \(B\). So if \(A\) has dimension \(m \times n\) then \(B\) needs to have dimension \(n \times r\) since matrix multiplication implies that we form the ‘inner product’ of each row of \(A\) with each column of \(B\). This results in a new matrix of dimension \(m \times r\). Here is an example. Given matrices \(A = \left[ \begin{array}{ccc} 12 & 3 &6 \\ 9 &-1 & -4 \end{array} \right] \text{ and } B = \left[ \begin{array}{cc} 7 & 8 \\ -2 &0 \\ 1 & 11 \end{array} \right]\)
the product of \(A\) \((2 \times 3)\) and \(B\) of dimension \((3 \times 2)\) is a matrix \(C\) with dimension \((2 \times 2)\):
\[AB = \left[ \begin{array}{cc} 12*7+3*(-2)+6*1 & 12*8+3*0+6*11 \\ 9*7+(-1)*(-2)+(-4)*1 & 9*8+(-1)*0+(-4)*11 \end{array} \right].\]
In Python matrix multiplication is achieved using the command dot() or the @ operator. This last operator is new and will only work with Python 3.5 and newer.
A = np.array([[12,3,6],[3,-1,-4]])
B = np.array([[7,8],[-2,0],[1,11]])
C = np.dot(A,B) # matrix multiplication
D = A @ B # new matrix multiplication operator
print("A=", A)
print("B=", B)
print("C=", C)
print("D=", D)A= [[12 3 6]
[ 3 -1 -4]]
B= [[ 7 8]
[-2 0]
[ 1 11]]
C= [[ 84 162]
[ 19 -20]]
D= [[ 84 162]
[ 19 -20]]Vectors are simply \(n \times 1\) or \(1 \times n\) dimensional matrices so that the same rules as above apply for multiplying a matrix with a vector.
A = np.array([[12,3,6],[3,-1,-4]])
b = np.array([7,-2,1])
C = np.dot(A,b) # matrix multiplied by vector: (m x n) x (n x 1) = (m x 1)
D = A @ b # new matrix multiplication operator
print("A=", A)
print("b=", b)
print("C=", C)
print("D=", D)A= [[12 3 6]
[ 3 -1 -4]]
b= [ 7 -2 1]
C= [84 19]
D= [84 19]