【机器学习-02】矩阵基础运算---numpy操作

  在机器学习-01中，我们介绍了关于机器学习的一般建模流程，并且在基本没有数学公式和代码的情况下，简单介绍了关于线性回归的一般实现形式。不过这只是在初学阶段、为了不增加基础概念理解难度所采取的方法，但所有的技术最终都是为了解决实际问题的，因此，接下来，我们就在之前的基础上更进一步，从一个更加严谨的理论体系出发、来尝试进行一种更加贴合实际应用所采用的一般方法的建模方法的学习。

import numpy as np
import pandas as pd

一、NumPy矩阵运算基础

  在进入到本节正式内容之前，我们需要先补充一些矩阵相关基础概念，以及矩阵运算的基本方法。
  在机器学习基础阶段，需要掌握的矩阵及线性代数基本理论包括：

• 矩阵的形变及特殊矩阵的构造方法：包括矩阵的转置、对角矩阵的创建、单位矩阵的创建、上/下三角矩阵的创建等；
• 矩阵的基本运算：包括矩阵乘法、向量内积、矩阵和向量的乘法等；
• 矩阵的线性代数运算：包括矩阵的迹、矩阵的秩、逆矩阵的求解、伴随矩阵和广义逆矩阵等；
• 矩阵分解运算：特征分解、奇异值分解和SVD分解等。
  本节将先介绍前三部分内容，矩阵分解部分内容将在后续补充。

1.NumPy中的矩阵表示

  在NumPy中，二维数组(array)和matrix类型对象都可以用于表示矩阵，并且也都具备矩阵的代数学方法。

• 利用数组创建矩阵

A = np.array([[1, 2], [1, 1]])Aarray([[1, 2],[1, 1]])type(A)numpy.ndarray

• 利用mat创建矩阵

AM = np.mat(A)
AMmatrix([[1, 2],[1, 1]])type(AM)
numpy.matrix

关于两种对象类型的选取，此处进行简单说明：

• NumPy中的matrix类型对象和MATLAB中的matrix类型等价，和NumPy中数组类型对象底层基本结构不同；
• 在NumPy中，针对大规模数据，数组类型对象的计算速度要快于矩阵类型对象；
• 矩阵类型对象可以通过运算符直接进行矩阵乘法，而二维数组要进行矩阵乘法（及其他矩阵运算），则必须要使用包括linalg（线性代数运算）模块在内的相关函数。

AM * AMmatrix([[3, 4],[2, 3]])A.dot(A)
array([[3, 4],[2, 3]])# 新版NumPy也支持使用符号进行矩阵乘法
A @ A
array([[3, 4],[2, 3]])

为了执行更高效的计算、以及确保代码整体基本对象类型统一，课程如无说明，将统一使用二维数组表示矩阵。

2.NumPy中特殊矩阵构造方法

  在实际线性代数运算过程中，经常涉及一些特殊矩阵，如单位矩阵、对角矩阵等，相关创建方法如下：

函数	描述
a.T	数组a转置
np.eye(n)	创建包含n个分量的单位矩阵
np.diag(a1)	以a1中各元素，创建对角矩阵
np.triu(a)	取矩阵a中的上三角矩阵
np.tril(a)	取矩阵a中的下三角矩阵

下面将展示一些矩阵的基本操作：

# 创建一个2*3的矩阵
a1 = np.arange(1, 7).reshape(2, 3)a1
array([[1, 2, 3],[4, 5, 6]])# 转置
a1.Tarray([[1, 4],[2, 5],[3, 6]])

矩阵的转置就是每个元素行列位置互换

# 创建单位矩阵
np.eye(3)
array([[1., 0., 0.],[0., 1., 0.],[0., 0., 1.]])

单位矩阵之所以被称为“单位”，核心原因在于单位矩阵和任何矩阵相乘，都将返回原矩阵。

a = np.arange(5)a
array([0, 1, 2, 3, 4])np.diag(a)
array([[0, 0, 0, 0, 0],[0, 1, 0, 0, 0],[0, 0, 2, 0, 0],[0, 0, 0, 3, 0],[0, 0, 0, 0, 4]])# 对角线向上偏移一位
np.diag(a, 1)
array([[0, 0, 0, 0, 0, 0],[0, 0, 1, 0, 0, 0],[0, 0, 0, 2, 0, 0],[0, 0, 0, 0, 3, 0],[0, 0, 0, 0, 0, 4],[0, 0, 0, 0, 0, 0]])# 对角线向下偏移一位
np.diag(a, -1)
array([[0, 0, 0, 0, 0, 0],[0, 0, 0, 0, 0, 0],[0, 1, 0, 0, 0, 0],[0, 0, 2, 0, 0, 0],[0, 0, 0, 3, 0, 0],[0, 0, 0, 0, 4, 0]])a1 = np.arange(9).reshape(3, 3)
a1
array([[0, 1, 2],[3, 4, 5],[6, 7, 8]])# 取上三角矩阵
np.triu(a1)
array([[0, 1, 2],[0, 4, 5],[0, 0, 8]])# 上三角矩阵向左下偏移一位
np.triu(a1, -1)
array([[0, 1, 2],[3, 4, 5],[0, 7, 8]])# 上三角矩阵向右上偏移一位
np.triu(a1, 1)
array([[0, 1, 2],[0, 0, 5],[0, 0, 0]])# 下三角矩阵
np.tril(a1)
array([[0, 0, 0],[3, 4, 0],[6, 7, 8]])

3.NumPy中矩阵基本运算

  由于NumPy中我们使用二维数组来表述矩阵，因此二维数组也就具备了数组和矩阵的两重属性。其中数组属性决定的基本运算相对简单，基础运算（如加减乘除）就是对应位置元素进行逐元素计算，而矩阵属性决定的运算则稍显复杂，当然矩阵的相关线性代数运算将在下一小节讨论，在基础运算上，矩阵和数组核心的区别在于乘法运算。
  当然，从另一个角度考虑，其实对于向量和矩阵这种具备一定结构的对象，有很多种容易混淆的计算规则。对于常用的计算规则，我们通过将其划分成三类以帮助大家理解：

描述	解释/函数
逐元素相乘	向量、矩阵通用
每个对应位置元素相乘	*
逐元素相乘后相加	也被称为点积（内积），向量，矩阵通用
向量点积	vdot、dot、inner
矩阵点积	vdot
矩阵乘法	代数学意义的矩阵相乘
矩阵乘法	dot、matmul、@

• * :逐元素相乘

a = np.arange(4)
a
array([0, 1, 2, 3])a * a
array([0, 1, 4, 9])A = a.reshape(2, 2)
A
array([[0, 1],[2, 3]])A * A
array([[0, 1],[4, 9]])

• 向量点积
    所谓点积（也被称为内积），指的是向量或矩阵对应位置元素相乘后相加。向量点积有三种实现方法，分别是dot、vdot和ineer。

np.dot(a, a)
14a.dot(a)
14(a * a).sum()
14np.vdot(a, a)
14np.inner(a, a)
14

• 矩阵点积
    值得注意的是，矩阵内积只有vdot一种方式实现。

A
array([[0, 1],[2, 3]])np.vdot(A, A)
14(A * A).sum()
14

注意，高维数组的inner并不是内积，而是一种类似tensordot的沿着尾轴实现和积的计算过程，该方法并不通用，此处暂不做介绍。

• 矩阵乘法
    NumPy中，我们可以使用诸多方法实现矩阵乘法，包括dot、@、matmul等。

a1 = np.arange(1, 7).reshape(2, 3)
a1
array([[1, 2, 3],[4, 5, 6]])a2 = np.arange(1, 10).reshape(3, 3)
a2
array([[1, 2, 3],[4, 5, 6],[7, 8, 9]])# 矩阵乘法
np.matmul(a1, a2)
array([[30, 36, 42],[66, 81, 96]])

此处也简单回顾矩阵乘法运算，上述相乘过程如下所示：

值得注意的是，矩阵相乘要求左乘矩阵列数和右乘矩阵行数相同，而内积计算过程则严格要求两个向量/矩阵形状完全一致。

4.NumPy中矩阵代数运算

  如果说矩阵的基本运算是矩阵基本性质，那么矩阵的线性代数运算，则是我们利用矩阵数据类型在求解实际问题过程中经常涉及到的线性代数方法，具体相关函数如下：

矩阵的线性代数运算

函数	描述
np.trace(A)	矩阵的迹
np.linalg.matrix_rank(A)	矩阵的秩
np.linalg…det(A)	计算矩阵A的行列式
np.linalg.inv(A)	矩阵求逆

同时，由于线性代数所涉及的数学基础知识较多，从实际应用的角度出发，我们将有所侧重的介绍实际应用过程中需要掌握的相关内容，并通过本节末尾的实际案例，来加深线性代数相关内容的理解。

NumPy中的linalg是linear algebra（线性代数）的简写，也是NumPy中保存线性代数相关计算函数的模块。

• 矩阵的迹（trace）

  矩阵的迹的运算相对简单，就是矩阵对角线元素之和，在NumPy中，可以使用trace函数进行计算。

A = np.array([[1, 2], [4, 5]])
A
array([[1, 2],[4, 5]])np.trace(A)
6

当然，对于矩阵的迹来说，计算过程不需要是方正

B = np.arange(1, 7).reshape(2, 3)
B
array([[1, 2, 3],[4, 5, 6]])np.trace(B)
6

• 矩阵的秩(rank)
    矩阵的秩（rank），是指矩阵中行或列的极大线性无关数，且矩阵中行、列极大无关数总是相同的，任何矩阵的秩都是唯一值，满秩指的是方阵（行数和列数相同的矩阵）中行数、列数和秩相同，满秩矩阵有线性唯一解等重要特性，而其他矩阵也能通过求解秩来降维，同时，秩也是奇异值分解等运算中涉及到的重要概念。

所谓线性相关，其实也就是线性表示，如果$y=wx+b$，我们则称y可以由x线性表示，二者线性相关，反之则线性无关。类似，如果$y=w_1{x}_1{w}_2{x}_2+b$，则我们称y可以由$x_1、x_2$线性表示，y与$x_1、x_2$线性相关。

• matrix_rank计算矩阵的秩

A = np.array([[1, 3, 4], [2, 1, 3], [1, 1, 2]])
A
array([[1, 3, 4],[2, 1, 3],[1, 1, 2]])np.linalg.matrix_rank(A)
2

对于矩阵A来说，第三列明显可以由第一列和第二列相加得出，因此极大线性无关组只有两列。

B = np.array([[1, 3, 4], [2, 1, 3], [1, 1, 10]])
B
array([[ 1,  3,  4],[ 2,  1,  3],[ 1,  1, 10]])np.linalg.matrix_rank(B)
3

• 矩阵的行列式(det)
    所谓行列式，我们可以简单将其理解为矩阵的一个基本性质或者属性，通过行列式的计算，我们能够知道矩阵是否可逆，从而可以进一步求解矩阵所对应的线性方程。当然，更加专业的解释，行列式的作为一个基本数学工具，实际上就是矩阵进行线性变换的伸缩因子。

对于任何一个n维方正，行列式计算过程如下：

更为简单的情况，如果对于一个2*2的矩阵，行列式的计算就是主对角线元素之积减去另外两个元素之积

A = np.array([[1, 2], [4, 5]])
A
array([[1, 2],[4, 5]])np.linalg.det(A)
-2.9999999999999996

A的秩计算过程如下：

对于行列式的计算，要求矩阵必须是方阵，也就是行列数必须一致。

B = np.arange(1, 7).reshape(2, 3)
B
array([[1, 2, 3],[4, 5, 6]])np.linalg.det(B)A = np.array([[1, 3, 4], [2, 1, 3], [1, 1, 2]])
A
array([[1, 3, 4],[2, 1, 3],[1, 1, 2]])np.linalg.det(A)
0.0

• 矩阵的逆
    对于满秩的方正来说，可以求其逆矩阵。从基本定义上来看，如果矩阵B和矩阵A相乘能够得到单位矩阵，即：
  $$B\cdot A=E$$

则称B为A的逆矩阵，也可将B写作$A^{-1}$。当然，逆矩阵的性质是相互的，我们也可称A为B的逆矩阵，或者A和B互为逆矩阵。

A = np.array([[1, 1], [3, 1]])
A
array([[1, 1],[3, 1]])

然后使用inverse函数进行逆矩阵求解

np.linalg.inv(A)
array([[-0.5,  0.5],[ 1.5, -0.5]])

简单试探逆矩阵的基本特性

A.dot(np.linalg.inv(A))
array([[1.00000000e+00, 1.11022302e-16],[2.22044605e-16, 1.00000000e+00]])

当然，对于逆矩阵，还有很多其他理解角度。例如，从方程组求解角度来看，逆矩阵的存在就代表着方程组存在唯一解，并且逆矩阵本身也是方程组求解的关键；从矩阵分解角度来看，逆矩阵是一种最为基础的矩阵分解的形式。关于这些相关内容，我们都将在后续课程中逐渐介绍。

另外，在本栏目内容中，我们还将介绍solve（方程组求解）、lstsq（最小二乘法）相关函数的使用。