71.numpy-applications-4.md

NumPy Applications - 4

Vectors

A vector is also called a directed quantity. It is a geometric object that has both magnitude and direction and satisfies the parallelogram law. The concept opposite to a vector is a scalar, which only has magnitude and in most cases has no direction. We usually use a line segment with an arrow to represent a vector. Vectors in the Cartesian plane are shown below. It should be noted that a vector expresses magnitude and direction, and it does not specify a start point or end point, so the same vector can be drawn at any position. For example, the two vectors $\small{\boldsymbol{w}}$ and $\small{\boldsymbol{u}}$ in the figure below are actually no different.

There are many algebraic ways to represent a vector. For a vector in two-dimensional space, the following notations are all acceptable:

$$ \boldsymbol{a} = \langle a_1, a_2 \rangle = (a_1, a_2) = \begin{pmatrix} a_1 \\ a_2 \end{pmatrix} = \begin{bmatrix} a_1 \\ a_2 \end{bmatrix} $$

The size of a vector is called its norm, and it is a scalar. For a vector in two-dimensional space, its norm can be calculated with the formula below:

$$ \lvert \boldsymbol{a} \rvert = \sqrt{a_{1}^{2} + a_{2}^{2}} $$

Notice that $\small{\lvert \boldsymbol{a} \rvert}$ here is not an absolute value. You can call it the 2-norm of vector $\small{\boldsymbol{a}}$. This is an example of symbol reuse in mathematics. The notation and concept above can be extended to $\small{n}$-dimensional space. We usually use $\small{\boldsymbol{R^{n}}}$ to represent $\small{n}$-dimensional space. The two-dimensional space just mentioned can be written as $\small{\boldsymbol{R^{2}}}$, and three-dimensional space can be written as $\small{\boldsymbol{R^{3}}}$. Although it is hard for us, who live in three-dimensional space, to imagine four-dimensional space or five-dimensional space, that does not stop us from discussing higher-dimensional space. In machine learning, we often call a training sample with $\small{n}$ features an $\small{n}$-dimensional vector.

Vector addition

Vectors of the same dimension can be added to produce a new vector. The operation is to add the corresponding components of the vectors, as shown below:

$$ \boldsymbol{u} = \begin{bmatrix} u_1 \\ u_2 \\ \vdots \\ u_n \end{bmatrix}, \quad \boldsymbol{v} = \begin{bmatrix} v_1 \\ v_2 \\ \vdots \\ v_n \end{bmatrix}, \quad \boldsymbol{u} + \boldsymbol{v} = \begin{bmatrix} u_1 + v_1 \\ u_2 + v_2 \\ \vdots \\ u_n + v_n \end{bmatrix} $$

Vector addition satisfies the parallelogram law. In other words, if vectors $\small{\boldsymbol{u}}$ and $\small{\boldsymbol{v}}$ form two adjacent sides of a parallelogram, their sum is the diagonal of that parallelogram, as shown below.

Scalar multiplication of a vector

A vector $\small{\boldsymbol{v}}$ can be multiplied by a scalar $\small{k}$. The operation is to multiply each component in the vector by that scalar, as shown below:

$$ \boldsymbol{v} = \begin{bmatrix} v_1 \\ v_2 \\ \vdots \\ v_n \end{bmatrix}, \quad k \cdot \boldsymbol{v} = \begin{bmatrix} k \cdot v_1 \\ k \cdot v_2 \\ \vdots \\ k \cdot v_n \end{bmatrix} $$

We can use NumPy arrays to represent vectors. Vector addition can be done by adding two arrays, and scalar multiplication of a vector can be done by multiplying an array by a scalar, so we will not repeat that here.

Dot product of vectors

The dot product is one of the most important operations between two vectors. The operation is to multiply corresponding components and then add the products together, so the result of a dot product is a scalar. Geometrically, it is equal to the product of the norms of the two vectors multiplied by the cosine of the angle between them.

$$ \boldsymbol{u} = \begin{bmatrix} u_1 \\ u_2 \\ \vdots \\ u_n \end{bmatrix}, \quad \boldsymbol{v} = \begin{bmatrix} v_1 \\ v_2 \\ \vdots \\ v_n \end{bmatrix} $$

$$ \boldsymbol{u} \cdot \boldsymbol{v} = \sum_{i=1}^{n}{u_iv_i} = \lvert \boldsymbol{u} \rvert \lvert \boldsymbol{v} \rvert cos\theta $$

Suppose we use a three-dimensional vector to represent a user's preference for comedy, romance, and action movies. We use numbers from 1 to 5 to represent the degree of liking, where 5 means like very much, 4 means like, 3 means neutral, 2 means dislike, and 1 means strongly dislike. Then the vector below represents a user who likes comedy very much, strongly dislikes romance, and feels neutral about action movies.

$$ \boldsymbol{u} = \begin{pmatrix} 5 \\ 1 \\ 3 \end{pmatrix} $$

Now suppose two movies are released. One is a romantic comedy, and the other is a comedy action movie. We also represent these two movies with 3-dimensional vectors, as shown below:

$$ \boldsymbol{m_1} = \begin{pmatrix} 4 \\ 5 \\ 1 \end{pmatrix}, \quad \boldsymbol{m_2} = \begin{pmatrix} 5 \\ 1 \\ 5 \end{pmatrix} $$

If we need to recommend one movie to this user, which one should we recommend? We can compute the dot product of the user vector $\small{\boldsymbol{u}}$ with the movie vectors $\small{\boldsymbol{m_{1}}}$ and $\small{\boldsymbol{m_{2}}}$, and then divide by the vector norms to get the cosine of the angle between the vectors. The closer the cosine value is to 1, the closer the angle is to 0 degrees, which means the two vectors are more similar. Obviously, we should recommend the movie whose vector is more similar to the user's viewing preference.

$$ cos\theta_1 = \frac{\boldsymbol{u} \cdot \boldsymbol{m1}}{|\boldsymbol{u}||\boldsymbol{m1}|} \approx \frac{4 \times 5 + 5 \times 1 + 3 \times 1}{5.92 \times 6.48} \approx 0.73 $$

$$ cos\theta_2 = \frac{\boldsymbol{u} \cdot \boldsymbol{m2}}{|\boldsymbol{u}||\boldsymbol{m2}|} \approx \frac{5 \times 5 + 1 \times 1 + 3 \times 5}{5.92 \times 7.14} \approx 0.97 $$

Some people may say that even by looking at it, the movie represented by vector $\small{\boldsymbol{m_{2}}}$ is obviously a better match for the user. That is true. For a three-dimensional vector, we can still rely on intuition to give the correct answer. But for an $\small{n}$-dimensional vector, when the value of $\small{n}$ becomes very large, can you still confidently rely on your eyes and intuition to give the right answer? The dot product of vectors can be calculated with the dot function, and the norm of a vector can be calculated with the norm function in NumPy's linalg module, as shown below.

u = np.array([5, 1, 3])
m1 = np.array([4, 5, 1])
m2 = np.array([5, 1, 5])
print(np.dot(u, m1) / (np.linalg.norm(u) * np.linalg.norm(m1)))  # 0.7302967433402214
print(np.dot(u, m2) / (np.linalg.norm(u) * np.linalg.norm(m2)))  # 0.9704311900788593

Cross product of vectors

In two-dimensional space, the cross product of two vectors is defined like this:

$$ \boldsymbol{A} = \begin{pmatrix} a_{1} \\ a_{2} \end{pmatrix}, \quad \boldsymbol{B} = \begin{pmatrix} b_{1} \\ b_{2} \end{pmatrix} $$

$$ \boldsymbol{A} \times \boldsymbol{B} = \begin{vmatrix} a_{1} \quad a_{2} \\ b_{1} \quad b_{2} \end{vmatrix} = a_{1}b_{2} - a_{2}b_{1} $$

For three-dimensional space, the cross product of two vectors is itself a vector:

$$ \boldsymbol{A} = \begin{pmatrix} a_{1} \\ a_{2} \\ a_{3} \end{pmatrix}, \quad \boldsymbol{B} = \begin{pmatrix} b_{1} \\ b_{2} \\ b_{3} \end{pmatrix} $$

$$ \boldsymbol{A} \times \boldsymbol{B} = \begin{vmatrix} \boldsymbol{\hat{i}} \quad \boldsymbol{\hat{j}} \quad \boldsymbol{\hat{k}} \\ a_{1} \quad a_{2} \quad a_{3} \\ b_{1} \quad b_{2} \quad b_{3} \end{vmatrix} $$

Since the result of a cross product is a vector, the results of $\small{\boldsymbol{A} \times \boldsymbol{B}}$ and $\small{\boldsymbol{B} \times \boldsymbol{A}}$ are not the same. In fact:

$$ \boldsymbol{A} \times \boldsymbol{B} = -(\boldsymbol{B} \times \boldsymbol{A}) $$

In NumPy, we can use the cross function to calculate the cross product of vectors, as shown below.

print(np.cross(u, m1))  # [-14   7  21]
print(np.cross(m1, u))  # [ 14  -7 -21]

Determinants

A determinant is usually written as $\small{det(\boldsymbol{A})}$ or $\small{\lvert \boldsymbol{A} \rvert}$, where $\small{\boldsymbol{A}}$ is an $\small{n}$-order square matrix. A determinant can be seen as the idea of signed area or signed volume extended to general Euclidean space, or you can say it describes the effect of a linear transformation on "volume". The idea of determinant first appeared when solving systems of linear equations. By the late seventeenth century, works by Seki Takakazu and Leibniz had already used determinants to decide the number and form of solutions of linear systems. From the eighteenth century on, determinants started to be studied as an independent mathematical concept, and after the nineteenth century, determinant theory became more complete.

Properties of determinants

Determinants come from vectors, so what determinants explain is really the properties of vectors.

1. If one row or one column of $\small{det(\boldsymbol{A})}$ is all 0, then $\small{det(\boldsymbol{A}) = 0}$.
2. If one row or one column of $\small{det(\boldsymbol{A})}$ has a common factor $\small{k}$, then $\small{k}$ can be taken out, and $\small{det(\boldsymbol{A}) = k \cdot det(\boldsymbol{A^{'}})}$.
3. If every element in one row or one column of $\small{det(\boldsymbol{A})}$ is the sum of two numbers, this determinant can be split into the sum of two determinants.
4. If two rows or two columns are proportional, then $\small{det(\boldsymbol{A}) = 0}$.
5. If swapping two rows or two columns gives $\small{det(\boldsymbol{A^{'}})}$, then $\small{det(\boldsymbol{A}) = -det(\boldsymbol{A^{'}})}$.
6. If we add $\small{k}$ times one row or one column to another row or column, the value of the determinant does not change.
7. If we swap rows and columns in a determinant, the value of the determinant does not change.
8. The determinant of the product of square matrices $\small{\boldsymbol{A}}$ and $\small{\boldsymbol{B}}$ is equal to the product of their determinants, that is $\small{det(\boldsymbol{A}\boldsymbol{B}) = det(\boldsymbol{A})det(\boldsymbol{B})}$. In a special case, if every row of a matrix is multiplied by the constant $\small{r}$, the determinant becomes $\small{r^n}$ times the original one.
9. If $\small{\boldsymbol{A}}$ is an invertible matrix, then $\small{det(\boldsymbol{A}^{-1}) = (det(\boldsymbol{A}))^{-1}}$.

Computing determinants

The general formula of an $\small{n}$-order determinant is shown below:

$$ det(\boldsymbol{A})=\sum_{n!} \pm {a_{1\alpha}a_{2\beta} \cdots a_{n\omega}} $$

For a second-order determinant, the formula above becomes:

$$ \begin{vmatrix} a_{11} \quad a_{12} \\ a_{21} \quad a_{22} \end{vmatrix} = a_{11}a_{22} - a_{12}a_{21} $$

For a third-order determinant, the formula above becomes:

$$ \begin{vmatrix} a_{11} \quad a_{12} \quad a_{13} \\ a_{21} \quad a_{22} \quad a_{23} \\ a_{31} \quad a_{32} \quad a_{33} \end{vmatrix} = a_{11}a_{22}a_{33} + a_{12}a_{23}a_{31} + a_{13}a_{21}a_{32} - a_{11}a_{23}a_{32} - a_{12}a_{21}a_{33} - a_{13}a_{22}a_{31} $$

Higher-order determinants can be expanded into several lower-order determinants by cofactors, as shown below:

$$ det(\boldsymbol{A})=a_{11}C_{11}+a_{12}C_{12}+ \cdots +a_{1n}C_{1n} = \sum_{i=1}^{n}{a_{1i}C_{1i}} $$

Here, $\small{C_{11}}$ is the determinant formed by the remaining part after removing the row and column where $\small{a_{11}}$ is. The others are understood in the same way.

Matrices

A matrix is a rectangular array made by arranging a series of elements. The elements in a matrix can be numbers, symbols, or mathematical formulas. Matrices can perform addition, subtraction, scalar multiplication, transpose, matrix multiplication, and other operations, as shown below.

One operation worth mentioning is matrix multiplication. It is defined only when the number of columns of the first matrix $\small{\boldsymbol{A}}$ is equal to the number of rows of the second matrix $\small{\boldsymbol{B}}$. If $\small{\boldsymbol{A}}$ is an $\small{m \times n}$ matrix and $\small{\boldsymbol{B}}$ is an $\small{n \times k}$ matrix, then their product is an $\small{m \times k}$ matrix.

For example:

$$ \begin{bmatrix} 1 & 0 & 2 \\ -1 & 3 & 1 \end{bmatrix} \times \begin{bmatrix} 3 & 1 \\ 2 & 1 \\ 1 & 0 \end{bmatrix} = \begin{bmatrix} 5 & 1 \\ 4 & 2 \end{bmatrix} $$

Matrix multiplication satisfies associativity and distributivity over matrix addition:

Associativity: $\small{(\boldsymbol{AB})\boldsymbol{C} = \boldsymbol{A}(\boldsymbol{BC})}$.

Left distributive law: $\small{(\boldsymbol{A} + \boldsymbol{B})\boldsymbol{C} = \boldsymbol{AC} + \boldsymbol{BC}}$.

Right distributive law: $\small{\boldsymbol{C}(\boldsymbol{A} + \boldsymbol{B}) = \boldsymbol{CA} + \boldsymbol{CB}}$.

Matrix multiplication does not satisfy the commutative law. In general, the product $\small{\boldsymbol{AB}}$ may exist while $\small{\boldsymbol{BA}}$ may not. Even if both exist, most of the time $\small{\boldsymbol{AB} \neq \boldsymbol{BA}}$.

One basic application of matrix multiplication is in systems of linear equations. A system of linear equations can be written in vector form as:

$$ \boldsymbol{Ax} = \boldsymbol{b} $$

where $\small{\boldsymbol{A}}$ is the coefficient matrix, $\small{\boldsymbol{x}}$ is the vector of unknowns, and $\small{\boldsymbol{b}}$ is the constant vector.

Matrices are also a convenient way to express linear transformations. If the discussion above feels hard to understand, I recommend the video series The Essence of Linear Algebra, which can give you a better understanding of linear algebra.

Matrix objects

NumPy provides a module for linear algebra and also a matrix type for representing matrices. Of course, we can also use a two-dimensional array to represent a matrix. Officially, using the matrix class is not recommended. It is better to use two-dimensional arrays, and the matrix class may be removed in future versions. Still, with these wrapped classes and functions, we can perform many matrix operations very easily.

We can create matrix objects with the following code:

m1 = np.matrix('1 2 3; 4 5 6')
m1

Note: The matrix constructor can accept either an array-like object or a string.

Output:

matrix([[1, 2, 3],
        [4, 5, 6]])

m2 = np.asmatrix(np.array([[1, 1], [2, 2], [3, 3]]))
m2

Note: The asmatrix function can also be replaced with the mat function. They are actually the same function.

Output:

matrix([[1, 1],
        [2, 2],
        [3, 3]])

m1 * m2

Output:

matrix([[14, 14],
        [32, 32]])

Note: Be careful about the difference between multiplication of matrix objects and ndarray objects. The * operator for a matrix object means matrix multiplication. If you want to do matrix multiplication with two two-dimensional arrays, you should use the @ operator or the matmul function, not the * operator.

Useful properties of matrix objects are shown below.

Property	Description
A	Get the corresponding ndarray object
A1	Get the flattened ndarray object
I	The inverse matrix of an invertible matrix
T	The transpose of the matrix
H	The conjugate transpose of the matrix
shape	The shape of the matrix
size	The number of elements in the matrix

The methods of matrix objects are almost the same as the methods of ndarray objects that we talked about earlier, so I will not repeat them here.

Linear algebra module

NumPy's linalg module contains a group of standard matrix decomposition operations, and functions such as inverse and determinant. The table below lists some commonly used linear-algebra-related functions in numpy and linalg.

Function	Description
diag	Return the diagonal elements of a square matrix as a one-dimensional array, or turn a one-dimensional array into a square matrix with zeros off the diagonal
matmul	Matrix multiplication
trace	Compute the sum of diagonal elements
norm	Compute the norm of a matrix or vector
det	Compute the determinant
matrix_rank	Compute the rank of a matrix
eig	Compute eigenvalues and eigenvectors
inv	Compute the inverse of a non-singular square matrix
pinv	Compute the Moore-Penrose generalized inverse
qr	QR decomposition
svd	Singular value decomposition
solve	Solve the linear system $\small{\boldsymbol{Ax}=\boldsymbol{b}}$ where $\small{\boldsymbol{A}}$ is a square matrix
lstsq	Compute the least-squares solution of $\small{\boldsymbol{Ax}=\boldsymbol{b}}$

Below we simply try several of the functions above. Let us first try finding an inverse matrix.

m3 = np.array([[1., 2.], [3., 4.]])
m4 = np.linalg.inv(m3)
m4

Output:

array([[-2. ,  1. ],
       [ 1.5, -0.5]])

np.around(m3 @ m4)

Note: The around function rounds array elements. By default, the number of decimal places is 0.

Output:

array([[1., 0.],
       [0., 1.]])

Note: A matrix multiplied by its inverse matrix gives the identity matrix.

Let us compute the determinant:

m5 = np.array([[1, 3, 5], [2, 4, 6], [4, 7, 9]])
np.linalg.det(m5)

Output:

2

Let us compute the rank of the matrix:

np.linalg.matrix_rank(m5)

Output:

3

Let us solve a linear system:

$$ \begin{cases} x_1 + 2x_2 + x_3 = 8 \\\ 3x_1 + 7x_2 + 2x_3 = 23 \\\ 2x_1 + 2x_2 + x_3 = 9 \end{cases} $$

For the system above, we can write it in matrix form as:

$$ \boldsymbol{A} = \begin{bmatrix} 1 & 2 & 1 \\\ 3 & 7 & 2 \\\ 2 & 2 & 1 \end{bmatrix}, \quad \boldsymbol{x} = \begin{bmatrix} x_1 \\\ x_2 \\\ x_3 \end{bmatrix}, \quad \boldsymbol{b} = \begin{bmatrix} 8 \\\ 23 \\\ 9 \end{bmatrix} $$

$$ \boldsymbol{Ax} = \boldsymbol{b} $$

The condition for a system of linear equations to have a unique solution is that the rank of the coefficient matrix $\small{\boldsymbol{A}}$ is equal to the rank of the augmented matrix $\small{\boldsymbol{Ab}}$, and both are equal to the number of unknowns.

A = np.array([[1, 2, 1], [3, 7, 2], [2, 2, 1]])
b = np.array([8, 23, 9]).reshape(-1, 1)
print(np.linalg.matrix_rank(A))
print(np.linalg.matrix_rank(np.hstack((A, b))))

Note: When reshaping an array object, if one argument is -1, the number of elements in that dimension is computed automatically from the total number of elements and the other dimension sizes.

Output:

3
3

np.linalg.solve(A, b)

Output:

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

Note: The result above shows that the solution to the linear system is $\small{x_1 = 1, x_2 = 2, x_3 = 3}$.

There is another way to solve the linear system. You can stop and think about why:

$$ \boldsymbol{x} = \boldsymbol{A}^{-1} \cdot \boldsymbol{b} $$

np.linalg.inv(A) @ b

Output:

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

Polynomials

Besides arrays, NumPy also wraps a data type for doing polynomial operations. A polynomial is a sum of products of coefficients and integer powers of a variable, in the form:

$$ f(x)=a_nx^n + a_{n-1}x^{n-1} + \cdots + a_1x^{1} + a_0x^{0} $$

Before NumPy 1.4, we could use the poly1d type to represent polynomials. It is still available now, but NumPy officially provides the newer numpy.polynomial module. In addition to ordinary power-series polynomials, it also supports Chebyshev polynomials, Laguerre polynomials, and others.

Creating polynomial objects

Create a poly1d object, for example $\small{f(x)=3x^{2}+2x+1}$:

p1 = np.poly1d([3, 2, 1])
p2 = np.poly1d([1, 2, 3])
print(p1)
print(p2)

Output:

   2
3 x + 2 x + 1
   2
1 x + 2 x + 3

Operations on polynomials

Get the coefficients of a polynomial

print(p1.coefficients)
print(p2.coeffs)

Output:

[3 2 1]
[1 2 3]

Four basic operations between two polynomials

print(p1 + p2)
print(p1 * p2)

Output:

   2
4 x + 4 x + 4
   4     3      2
3 x + 8 x + 14 x + 8 x + 3

Substitute a value of $\small{x}$ and evaluate the polynomial

print(p1(3))
print(p2(3))

Output:

34
18

Differentiate and integrate a polynomial

print(p1.deriv())
print(p1.integ())

Output:

6 x + 2
   3     2
1 x + 1 x + 1 x

Find the roots of a polynomial

For example, given $\small{f(x)=x^2+3x+2}$, the roots of the polynomial are the solutions of the quadratic equation $\small{x^2+3x+2=0}$.

p3 = np.poly1d([1, 3, 2])
print(p3.roots)

Output:

[-2. -1.]

If we use the Polynomial class from the numpy.polynomial module to represent polynomial objects, then the corresponding operations are as follows:

from numpy.polynomial import Polynomial

p3 = Polynomial((2, 3, 1))
print(p3)           # print the polynomial
print(p3(3))        # let x=3 and evaluate the polynomial
print(p3.roots())   # compute the roots of the polynomial
print(p3.degree())  # get the degree of the polynomial
print(p3.deriv())   # differentiate
print(p3.integ())   # compute the indefinite integral

Output:

2.0 + 3.0·x + 1.0·x²
20.0
[-2. -1.]
2
3.0 + 2.0·x
0.0 + 2.0·x + 1.5·x² + 0.33333333·x³

Least-squares solution

The Polynomial class also has a class method named fit, which can find the least-squares solution of a polynomial. The so-called least-squares solution is to use the least-squares method to find the coefficients of the function that best matches the data by minimizing the sum of squared errors. Suppose the polynomial is $\small{f(x)=ax+b}$, then the least-squares solution is the values of $\small{a}$ and $\small{b}$ that make the residual sum of squares $\small{RSS}$ below as small as possible.

$$ RSS = \sum_{i=0}^{k}{(f(x_i) - y_i)^{2}} $$

For example, suppose we want to use collected historical data of monthly income and online-shopping spending to build a prediction model, so that we can predict a person's online-shopping spending amount from that person's monthly income. The income and spending data we collected are stored in the two arrays below.

x = np.array([
    25000, 15850, 15500, 20500, 22000, 20010, 26050, 12500, 18500, 27300,
    15000,  8300, 23320,  5250,  5800,  9100,  4800, 16000, 28500, 32000,
    31300, 10800,  6750,  6020, 13300, 30020,  3200, 17300,  8835,  3500
])
y = np.array([
    2599, 1400, 1120, 2560, 1900, 1200, 2320,  800, 1650, 2200,
     980,  580, 1885,  600,  400,  800,  420, 1380, 1980, 3999,
    3800,  725,  520,  420, 1200, 4020,  350, 1500,  560,  500
])

We can first draw a scatter plot to see whether the two groups of data have a positive correlation or a negative correlation. Positive correlation means larger values in array x also match larger values in array y, while negative correlation means larger values in array x match smaller values in array y.

import matplotlib.pyplot as plt

plt.figure(dpi=120)
plt.scatter(x, y, color='blue')
plt.show()

Output:

If we want to study the correlation of the two groups of data quantitatively, we can compute the covariance or correlation coefficient. The corresponding NumPy functions are cov and corrcoef.

np.corrcoef(x, y)

Output:

array([[1.        , 0.92275889],
       [0.92275889, 1.        ]])

Note: The correlation coefficient is a value between -1 and 1. The closer it is to 1, the stronger the positive correlation is. The closer it is to -1, the stronger the negative correlation is. If it is close to 0, it means the two groups of data do not have an obvious correlation. The correlation coefficient between monthly income and online-shopping spending above is 0.92275889, which means they have a strong positive correlation.

From the work above, we confirm that there is a strong positive correlation between income and online-shopping spending, so we use these data to build a regression model and find a straight line that can fit these data points well. Here we can use the fit method mentioned above. The code is shown below.

from numpy.polynomial import Polynomial

Polynomial.fit(x, y, deg=1).convert().coef

Note: deg=1 means the highest power in the regression model is first degree, so the model has the form $\small{y=ax+b}$. If you want a model like $\small{y=ax^2+bx+c}$, you need to set deg=2, and so on.

Output:

array([-2.94883437e+02,  1.10333716e-01])

According to the output above, our regression equation should be $\small{y=0.110333716x-294.883437}$ . We draw this regression equation on the scatter plot we just showed. The red points are our predicted values, and the blue points are the historical data, which are the real values.

import matplotlib.pyplot as plt

plt.scatter(x, y, color='blue')
plt.scatter(x, 0.110333716 * x - 294.883437, color='red')
plt.plot(x, 0.110333716 * x - 294.883437, color='darkcyan')
plt.show()

Output:

If we do not use the fit method of the Polynomial type, we can also do the same thing with the polyfit function provided by NumPy. Interested readers can study that by themselves.

Note: Some images in this chapter come from Wikipedia.