下载掌阅APP，畅读海量书库
立即打开
NumPy数组也可以用作矩阵，这是数学和计算编程的基础。简单地说， 矩阵 就是一个二维数组。矩阵在许多应用中都是核心，如几何变换和联立方程，而在统计学等其他领域也作为有用的工具出现。与其他数组相比，矩阵本身只有在我们为其赋予矩阵算术运算时才是独特的。矩阵具有元素级的加法和减法运算，就像NumPy数组一样。还有一种称为标量乘法的第三种运算，就是将矩阵的每个元素乘以一个常数，它是一种不同的矩阵乘法概念。我们将在后面看到，矩阵乘法与其他乘法概念有着根本的不同。
矩阵最重要的属性之一是其形状，其定义与NumPy数组完全相同。通常具有 m 行和 n 列的矩阵被描述为 m × n 矩阵。行数与列数相同的矩阵称为方阵，这些矩阵在向量和矩阵理论中起着特殊的作用。
单位矩阵 （大小为 n ）是一个 n × n 矩阵，其中（ i , i ）位置上的元素为1，当 i ≠ j 时，（ i , j ）位置上的元素为零。下面的数组创建例程可以通过指定 n 的值来生成 n × n 单位矩阵：
顾名思义，单位矩阵是一个特殊的矩阵，具有如下性质：任何矩阵与单位矩阵相乘等于其本身、两个逆矩阵的乘积为单位矩阵等。
与单位矩阵相关的术语和量有很多。我们在这里只介绍两个性质，因为稍后会用到它们。这两个性质分别是矩阵的转置（transpose）和方阵的迹（trace）。其中矩阵转置表示行和列互换；迹是方阵主对角线的元素之和，主对角线由元素 a i,i 构成，即沿矩阵左上角到右下角的对角线上的元素。
NumPy数组可以通过在array对象上调用transpose方法轻松地进行转置。事实上，由于这是一种常见的操作，数组还有一个便捷的属性T，它也返回矩阵的转置。转置将矩阵（数组）的形状顺序颠倒，使得行变为列，列变为行。例如，如果我们有一个3×2矩阵（3行，2列），那么它的转置将是一个2×3矩阵，如下例所示：
注意
transpose函数实际上并不修改底层数组中的数据，而是改变数组的形状和内部标志，该标志表示存储数值的顺序是行连续（C风格）还是列连续（F风格）。这使得操作成本非常低。
与矩阵相关的另一个有用的量是迹。如前面的代码所示，方阵 A 的迹被定义为沿矩阵左上角到右下角主对角线上的元素之和。迹的公式如下所示：
NumPy数组有一个trace方法，可以返回矩阵的迹：
也可以使用np.trace函数获取迹，该函数并不是数组的绑定操作。
矩阵乘法是对两个矩阵执行的运算，它保留了两个矩阵的一些结构和特性。形式上，假设我们有如下 l × m 矩阵 A 和 m × n 矩阵 B 两个矩阵：
矩阵 A 和 B 的矩阵积 C 是一个 l × n 矩阵，其（ p , q ）项由以下方程给出：
注意，为了定义两个矩阵的乘法，第一个矩阵的列数必须与第二个矩阵的行数相匹配。通常，我们用 AB 表示矩阵 A 和 B 的乘积（如果根据定义它们能够相乘的话）。矩阵乘法是一种特殊的运算，它不像大多数其他算术运算那样满足交换律：即使 AB 和 BA 都可以计算，它们也不一定相等。在实践中，这意味着矩阵乘法的顺序很重要。这源于矩阵代数作为线性映射表示的起源，其中乘法对应于函数的组合。
Python有一个为矩阵乘法保留的运算符@，它是在Python 3.5中新添加的。NumPy数组应用该运算符能够实现矩阵乘法。请注意，@与数组的逐元素乘法运算符*有本质区别：
单位矩阵是矩阵乘法的中性元素。也就是说，如果 A 是任意 k × m 矩阵，而 I 是 m × m 单位矩阵，那么 AI = A ；同样，如果 B 是 m × k 矩阵，那么 IB = B 。使用NumPy数组可以很容易地检验特定示例：
可以看出，输出的结果矩阵等于原始矩阵。如果我们颠倒 A 和 I 的顺序并执行乘法 IA ，结果也是相同的。在1.5.3节中，我们将讨论矩阵的逆：若矩阵 B 与矩阵 A 相乘得到单位矩阵，则称矩阵 B 是 A 的逆矩阵。
方阵的行列式在很多应用中都很重要，因为它与求解矩阵的逆有着紧密的联系。如果矩阵的行数和列数相等，则它是一个方阵。特别地，具有非零行列式的矩阵有一个（唯一的）逆矩阵，这可以转化为某些方程组的唯一解。矩阵的行列式是通过递归定义的。假设我们有一个通用的2×2矩阵，如下所示：
这个通用矩阵 A 的行列式由以下公式定义：
对于一般的 n × n 矩阵，其中 n ＞2，我们通过递归的形式定义行列式。对于1≤ i ， j ≤ n ， A 的第 i 个0子矩阵 A i,j 是从矩阵 A 中删除第 i 行和第 j 列后的结果。子矩阵 A i,j 是（ n -1）×（ n -1）矩阵，因此我们可以计算它的行列式。然后，我们将以下量定义为 A 的行列式：
实际上，上述方程中出现的索引1可以替换为任意1≤ i ≤ n ，结果是相同的。
在NumPy中，用于计算行列式的例程位于名为linalg的独立模块中。该模块包含许多线性代数（涵盖向量和矩阵代数的数学分支）的常用函数。用于计算方阵行列式的具体例程为det：
请注意，在计算行列式时出现了浮点舍入误差。
如果安装了SciPy包，它也提供了一个linalg模块，该模块扩展了NumPy的linalg。SciPy版本的linalg模块不仅包括额外的例程，而且它在编译时始终支持 基础线性代数子程序包 （BLAS）和 线性代数包 （LAPACK），而它们在NumPy版本中则是可选的。因此，如果更注重速度，根据NumPy的编译方式，使用SciPy版本可能更好。
n × n 矩阵 A 的 逆矩阵 是另一个（必然唯一的） n × n 矩阵 B ，使得 AB = BA = I ，其中 I 表示 n × n 单位矩阵，这里的乘法是矩阵乘法。并不是每个方阵都有逆矩阵，那些没有逆矩阵的矩阵有时被称为 奇异矩阵 。事实上，当且仅当一个矩阵的行列式不为0时，该矩阵才是非奇异的（即具有逆矩阵）。当 A 有逆矩阵时，它的逆矩阵通常用 A -1 表示。
用linalg模块中的inv例程可以计算矩阵的逆（如果存在）：
将矩阵 A 乘以其逆矩阵（在任一侧）得到的结果是2×2单位矩阵，依据这一事实我们可以验证inv例程给出的矩阵确实是 A 的逆矩阵：
由于计算逆矩阵的方式不同，在这些计算中会出现浮点误差，这些误差已经隐藏在以“Approximately”（约等于）开头的注释后面。
linalg包还包含许多其他方法，如用于计算矩阵各种范数的norm方法。它还包括以各种方式分解矩阵和求解方程组的函数。
除此之外，它还包含用于矩阵运算的指数函数expm、对数函数logm、正弦函数sinm、余弦函数cosm和正切函数tanm。请注意，这些函数与NumPy基本包中的标准exp、log、sin、cos和tan函数不同，后者是按逐元素方式执行相应函数的。相反，矩阵指数函数是使用矩阵的幂级数定义的：
这是对任意 n × n 矩阵 A 定义的，其中 A k 表示 A 的 k 次矩阵幂，即 A 矩阵连乘自己 k 次。请注意，这个“幂级数”总是在适当的意义上收敛。按照惯例，我们取 A 0 = I ，其中 I 是 n × n 单位矩阵。这完全类似于实数或复数的指数函数的常规幂级数定义，但用矩阵和矩阵乘法代替了数字和常规乘法。其他函数以类似的方式定义，但我们将略过这些详细信息。
在1.5.4节中，我们将看到矩阵及其理论在求解方程组领域的应用。
求解线性方程组是数学中研究矩阵的主要动机之一，这类问题在各种应用中频繁出现。我们从如下线性方程组开始：
这里， n 至少为2， a i,j 和 b i 是已知量， x i 是我们希望求解的未知数。
在求解这样的方程组之前，我们需要将问题转化为矩阵方程。这是通过将方程中的系数 a i,j 收集到一个 n × n 矩阵中，并利用矩阵乘法的性质将这个矩阵与方程组关联起来来实现的。因此，我们构造以下矩阵，矩阵元素为方程的系数：
然后，如果我们将 x 定义为包含 x i 值的未知（列）向量，将 b 定义为包含 b i 值的已知（列）向量，那么我们可以将方程组重写为以下单矩阵方程：
我们现在可以使用矩阵技术来求解这个矩阵方程。在这种情况下，我们将列向量视为 n ×1矩阵，因此上述方程中的乘法是矩阵乘法。为了求解该矩阵方程，我们使用linalg模块中的solve例程。为了说明该技术，我们将以求解以下方程组为例：
这个方程组有三个未知量： x 1 、 x 2 和 x 3 。首先，我们创建系数矩阵和向量 b 。由于我们使用NumPy来处理矩阵和向量，我们为矩阵 A 创建一个二维NumPy数组，并为 b 创建一个一维数组：
现在，可以使用linalg.solve例程求出方程组的解：
这确实是方程组的解，可以通过计算A@x并将结果与b数组进行比较来轻松验证。在这个计算中可能存在浮点舍入误差。
solve例程需要两个输入，即系数矩阵 A 和等式右侧的向量 b 。它使用的例程将矩阵 A 分解为更简单的矩阵，从而快速将问题简化为可通过简单替换求解的形式。这种求解矩阵方程的技术非常强大和高效，而且不太容易出现困扰其他方法的浮点舍入误差。例如，如果矩阵的逆是已知的，那么可以通过左乘 A 的逆矩阵来计算方程组的解。然而，这通常不如使用solve例程效果好，因为它可能更慢或者会导致更大的数值误差。
在我们使用的例子中，系数矩阵 A 是方阵。也就是说，方程的数量与未知数的数量相同。在这种情况下，当且仅当矩阵 A 的行列式不为0时，方程组才有唯一解。当 A 的行列式为0时，可能发生两种情况：一是方程组没有解，这时我们说方程组是不相容的；二是方程组可能有无穷多个解。方程组是相容还是不相容通常由向量 b 决定。例如，考虑以下方程组：
左边的方程组是相容的，并且有无穷多个解：例如，取 x =1和 y =1，或者取 x =0和 y =2都是解。右边的方程组是不相容的，没有解。在这两个方程组中，solve例程都会失败，因为系数矩阵是奇异的。
求解方程组时，不一定要求系数矩阵为方阵——例如，如果方程的数量多于未知数的数量（系数矩阵的行数多于列数）。这样的方程组称为“超定方程组”，只要它是相容的，就会有解。如果方程的数量少于未知数的数量，则称该方程组为“欠定方程组”。欠定方程组没有足够的信息来唯一确定所有的未知数，因此如果它是相容的，它通常有无穷多个解。遗憾的是，对于系数矩阵不是方阵的方程组，即使方程组确实有解，solve例程也无法找到解。
在1.5.5节中，我们将讨论特征值和特征向量，与你在前面看到的类似，它们是通过观察一种非常特殊的矩阵方程产生的。
考虑矩阵方程 Ax = λ x ，其中 A 是一个 n × n 方阵， x 是一个向量， λ 是一个系数。若存在 x 是此方程的解，则称 λ 为 A 的特征值，对应的向量 x 称为特征向量。特征值和对应的特征向量对矩阵 A 的信息进行编码，因此在许多涉及矩阵的应用中非常重要。
我们使用以下矩阵演示如何计算特征值和特征向量：
我们首先需要将其定义为NumPy数组：
linalg模块中的eig例程可以用来找到方阵的特征值和特征向量。该函数返回一个二元组（v,B），其中v是包含特征值的一维数组，B是一个二维数组，其列是相应的特征向量：
只有实数元素的矩阵完全可能具有复特征值和特征向量。因此，eig例程的返回类型有时会是复数类型，如complex32或complex64。在某些应用中，复特征值具有特殊含义，而在其他情况下，我们只考虑实特征值。
我们可以使用以下序列从eig例程的输出结果中提取特征值和特征向量：
eig例程返回的特征向量已经被归一化，使得它们的范数（长度）为1。（欧氏范数被定义为数组元素平方和的平方根。）我们可以通过使用linalg模块中的norm例程计算特征向量的范数来验证这一点：
最后，我们可以通过计算A @ x0，并检查它是否在浮点精度范围内等于lamb-da0*x0，来验证这些值确实满足特征值和特征向量的定义：
这里计算得到的范数表示方程 Ax = λ x 的左侧（lhs）和右侧（rhs）之间的距离。由于这个距离非常小（小数点后14位为零），我们可以相当确信它们实际上是相同的。事实上，这种距离不为零可能是由于浮点精度误差引起的。
求特征值和特征向量的理论过程是先通过解以下方程来求特征值 λ ：
这里， I 是适当的单位矩阵。左侧确定的方程是关于 λ 的多项式，称为矩阵 A 的特征多项式。然后，可以通过求解以下矩阵方程来求出与特征值 λ i 相对应的特征向量：
在实践中，这个过程效率有些低，还有其他策略可以更高效地通过数值计算求得特征值和特征向量。
我们只能计算方阵的特征值和特征向量，对于不是方阵的矩阵，这个定义没有意义。有一种称为奇异值的概念可以将特征值和特征向量推广到非方阵。为了做到这一点，我们必须做出的权衡是，我们必须计算两个向量 u 和 v 以及奇异值 σ ，然后求解以下方程：
如果 A 是 m × n 矩阵，那么 u 将有 n 个元素， v 将有 m 个元素。有趣的是， u 向量实际上是对称矩阵 A T A 的（正交归一化的）特征向量，其特征值为 σ 2 。根据这些值，我们可以利用之前的定义方程求出 v 向量。这将产生所有有趣的组合，但还有其他的向量 u 和 v ，使得 Au =0和 A T v =0。
奇异值及对应向量的效用来自 奇异值分解 （SVD），它将矩阵 A 写成以下乘积形式：
这里， U 为正交列， V 为正交行， Σ 为对角矩阵，通常的写法是数值沿主对角线递减。我们可以用稍微不同的方式写出这个公式，如下所示：
这就是说，任何矩阵都可以分解成外积的加权和——假设 u 和 v 是 n 行1列的矩阵，然后将 u 与向量 v 的转置进行矩阵乘法。
一旦完成了这种分解，我们就可以寻找特别小的 σ 值，这些值对矩阵贡献很小。如果我们放弃具有小 σ 值的项，那么我们就可以用更简单的表示法来有效地近似原始矩阵。这种技术在 主成分分析 （PCA）中得到了应用——例如，将复杂的高维数据集简化为对数据整体特征贡献最大的几个分量。
在Python中，我们可以使用linalg.svd函数来计算矩阵的SVD。该函数的工作方式与之前描述的eig例程类似，只是它返回的是分解后的三个分量：
该函数返回的数组形状分别为（2, 2）、（2,）和（4, 4）。顾名思义，U矩阵和VT矩阵是分解中出现的矩阵，而s是包含非零奇异值的一维向量。我们可以通过重构 Σ 矩阵并计算三个矩阵的乘积来检查分解是否正确：
请注意，除了矩阵的第一个元素外，矩阵几乎完全被重构了。左上角元素的值非常接近零——在浮点误差范围内可以被视为零。
我们构建矩阵 Σ 的方法相当不方便。SciPy版本的linalg模块包含一个特殊例程linalg.diagsvd，用于从一维奇异值数组重构该矩阵。该函数获取奇异值数组s和原始矩阵的形状，并构建具有适当形状的矩阵 Σ ：
（回想一下，SciPy软件包是以别名sp导入的。）现在，让我们改变一下方向，看看如何更有效地描述大多数元素都为零的矩阵，这就是所谓的稀疏矩阵。
前面所讨论的线性方程组，在整个数学领域，特别是在数学计算中是非常常见的。在许多应用中，系数矩阵可能会非常庞大，具有成千上万的行和列，并且很可能来自其他来源而不是简单地手动输入。在许多情况下，这也可能是一个稀疏矩阵，其中大多数元素为0。
如果一个矩阵有大量元素为零，则这个矩阵称为 稀疏矩阵 。矩阵有多少个元素为零才能称为稀疏矩阵，这并没有确切的定义。稀疏矩阵可以通过更高效的表示方式来存储，例如，只存储非零元素的索引（ i , j ）和相应的值 a i,j 。对于稀疏矩阵，有一整套专门的算法集合，可以在矩阵足够稀疏的情况下显著提高性能。
稀疏矩阵出现在许多应用中，并且通常遵循某种模式。特别是，求解 偏微分方程 （PDE）的几种技术都涉及求解稀疏矩阵方程（参见第3章），与网络相关的矩阵通常也是稀疏的。在sparse.csgraph模块中包含了与网络（图）相关的稀疏矩阵的附加例程。我们将在第5章中进一步讨论这些内容。
sparse模块包含几个不同的类，表示存储稀疏矩阵的不同方法。存储稀疏矩阵最基本的方法是存储三个数组，其中两个数组包含表示非零元素索引的整数，第三个数组包含相应的元素数据，这是coo_matrix类的格式。此外，还有 压缩稀疏列 格式（CSC，使用csc_matrix实现）和 压缩稀疏行 格式（CSR，使用csr_matrix实现），分别提供高效的列切片或行切片。在sparse模块中还有三个额外的稀疏矩阵类，包括dia_matrix类，它可以高效地存储非零元素沿对角线带出现的矩阵。
SciPy的sparse模块包含用于创建和处理稀疏矩阵的例程。我们使用以下导入语句从SciPy导入sparse模块：
稀疏矩阵可以由满矩阵（稠密矩阵）或其他类型的数据结构创建。这是使用特定格式（你希望使用这样的格式存储稀疏矩阵）的构造函数来完成的。
例如，我们可以使用以下命令将稠密矩阵存储为CSR格式：
如果手动生成稀疏矩阵，则该矩阵可能遵循某种模式，例如以下三对角矩阵：
这里，非零元素出现在对角线上以及对角线的两侧，并且每行的非零元素遵循相同的模式。要创建这样的矩阵，我们可以使用sparse中的一个数组创建例程，例如diags，这是一个创建具有对角线模式的矩阵的便捷例程：
这将创建一个如前所述的矩阵 T ，并将其存储为CSR格式的稀疏矩阵。第一个参数指定应出现在输出矩阵中的值，第二个参数是放置数值时相对于对角线的位置。因此，元组中的索引0表示对角线元素，-1表示在行中的对角线左侧，+1表示在行中的对角线右侧。shape关键字参数给出了所生成矩阵的维度，format指定了矩阵的存储格式。如果没有使用可选参数来设置格式，则将使用合理的默认值。可以使用toarray方法将数组T扩展为满矩阵（稠密矩阵）：
当矩阵很小时（如这里所示），稀疏矩阵求解函数与常规的求解函数之间的性能几乎没有差异。
一旦将矩阵以稀疏格式存储，我们就可以使用linalg中的子模块sparse求解例程。例如，我们可以使用此模块中的spsolve例程来求解矩阵方程。如果矩阵不是以稀疏格式提供的，spsolve例程会将矩阵转换为CSR或CSC格式，这可能会增加计算时间：
sparse.linalg模块还包含许多可以在NumPy（或SciPy）的linalg模块中找到的函数，这些函数接受稀疏矩阵而不是完整的NumPy数组，例如eig和inv。
至此，我们对Python及其生态系统中可用的基本数学工具的介绍就结束了。让我们总结一下我们所学到的内容。