科学计算基础软件包NumPy入门讲座：线性代数子模块

Posted 2022-05-16 天元浪子

文章目录

• 1. 计算矩阵的行列式
• 2. 求解逆矩阵
• 3. 计算特征向量和特征值
• 4. 矩阵的奇异值分解
• 5. 求解线性方程组

NumPy 的线性代数子模块（linalg）提供了 20 余个函数，用于求解行列式、逆矩阵、特征值、特征向量，以及矩阵分解等。SciPy 的线性代数子模块（同样名为 linalg）更为庞大，提供了超过一百个函数。两个 linalg 子模块的同名函数基本保持了相同的功能，有些函数可能略有差异。为了尽可能同时给出两个模块同名函数的应用示例，本节的代码同时导入两个子模块，一个命名为 sla，另一个命名为 nla。

>>> import numpy as np
>>> from scipy import linalg as sla
>>> from numpy import linalg as nla

1. 计算矩阵的行列式

行列式在本质上可以视为线性变换的伸缩因子，因此行列式是一个标量。如果一个方阵（行数和列数相等的矩阵）的行列式等于 0，则该方阵为奇异矩阵，否则为非奇异矩阵。计算矩阵的行列式虽然简单，但手工计算很容易出错，而使用 linalg.det( ) 函数来计算则是易如反掌。

>>> m = np.mat('0 1 2; 1 0 3; 4 -3 8')
>>> sla.det(m) # scipy.linalg
-1.9999999999999982
>>> nla.det(m) # numpy.linalg
-2.0

2. 求解逆矩阵

矩阵可逆的条件是非奇异，也就是行列式不等于 0。从数学的角度看，矩阵可逆或非奇异是好的属性，类似函数的可微、可导。尽管 matrix 对象本身有逆矩阵的属性，但用 linalg 子模块求解矩阵的逆，也是非常简单的。

>>> m = np.mat('0 1 2; 1 0 3; 4 -3 8')
>>> m.I # matrix对象的逆矩阵属性
matrix([[-4.5, 7. , -1.5],
 [-2. , 4. , -1. ],
 [ 1.5, -2. , 0.5]])
>>> m*m.I # 矩阵和其逆矩阵的乘积为单位矩阵
matrix([[1., 0., 0.],
 [0., 1., 0.],
 [0., 0., 1.]])
>>> sla.inv(m) # scipy.linalg
array([[-4.5, 7. , -1.5],
 [-2. , 4. , -1. ],
 [ 1.5, -2. , 0.5]])
>>> nla.inv(m) # numpy.linalg
matrix([[-4.5, 7. , -1.5],
 [-2. , 4. , -1. ],
 [ 1.5, -2. , 0.5]])

3. 计算特征向量和特征值

对于 n 阶矩阵 A，若存在标量 λ 和 n 维非零列向量 x，使得 Ax=λx，那么标量 λ 就称为矩阵A 的特征值，向量 x 就称为矩阵 A 的特征向量。n 阶矩阵存在 n 个特征值，每个特征值对应一个特征向量。

>>> A = np.mat('0 1 2; 1 0 3; 4 -3 8') # 生成3阶矩阵
>>> sla.eigvals(A) # 返回3个特征值
array([ 7.96850246+0.j, -0.48548592+0.j, 0.51698346+0.j])
>>> sla.eig(A) # 返回3个特征值和3个特征向量组成的元组
(array([ 7.96850246+0.j, -0.48548592+0.j, 0.51698346+0.j]), 
array([[ 0.26955165, 0.90772191, -0.74373492],
 [ 0.36874217, 0.24316331, -0.65468206],
 [ 0.88959042, -0.34192476, 0.13509171]]))
>>> nla.eigvals(A) # 返回3个特征值
array([ 7.96850246, -0.48548592, 0.51698346])
>>> nla.eig(A) # 返回3个特征值和3个特征向量组成的元组
(array([ 7.96850246, -0.48548592, 0.51698346]), 
matrix([[ 0.26955165, 0.90772191, -0.74373492],
 [ 0.36874217, 0.24316331, -0.65468206],
 [ 0.88959042, -0.34192476, 0.13509171]]))

4. 矩阵的奇异值分解

特征向量和特征值是矩阵最重要的特征，特征向量表示特征是什么，特征值表示这个特征有多重要。特征向量和特征值是通过对矩阵的特征分解获得的，不过特征分解只适用于方阵。实际应用中，很多矩阵并不是方阵，要想获取矩阵特征就要使用矩阵的奇异值分解。

对 m×n 阶矩阵进行奇异值分解，返回一个三元组：以左奇异向量为列的矩阵 U，形状为(m, m) 或 (m, k) ；按降序排列的奇异值向量 s，形状为 (k, )，其中 k=min(m, n) ；以右奇异向量为行的矩阵 V，形状为 (n, n) 或 (k, n)。

>>> A = np.mat(np.random.randint(0,10,(3,4)))
>>> A
matrix([[6, 7, 8, 1],
 [7, 4, 9, 6],
 [4, 6, 2, 1]])
>>> U, s, V = sla.svd(A)
>>> U.shape, s.shape, V.shape
((3, 3), (3,), (4, 4))
>>> U
array([[-0.63408037, -0.38759249, -0.66911445],
 [-0.68934803, 0.67538003, 0.26203265],
 [-0.35034465, -0.62740249, 0.69543134]])
>>> s
array([18.88807714, 5.14426409, 2.40355755])
>>> V
array([[-0.53109149, -0.49226941, -0.63412832, -0.27109764],
 [-0.02089797, -0.73402962, 0.33491192, 0.59042171],
 [ 0.2501572 , 0.22338451, -0.66724387, 0.66506116],
 [-0.80927528, 0.41106046, 0.20124835, 0.36824166]])
>>> U, s, V = nla.svd(A)
>>> U
matrix([[-0.63408037, -0.38759249, -0.66911445],
 [-0.68934803, 0.67538003, 0.26203265],
 [-0.35034465, -0.62740249, 0.69543134]])
>>> s
array([18.88807714, 5.14426409, 2.40355755])
>>> V
matrix([[-0.53109149, -0.49226941, -0.63412832, -0.27109764],
 [-0.02089797, -0.73402962, 0.33491192, 0.59042171],
 [ 0.2501572 , 0.22338451, -0.66724387, 0.66506116],
 [-0.80927528, 0.41106046, 0.20124835, 0.36824166]])

5. 求解线性方程组

求解线性方程组对中学生来说是一件轻松的事情，但是用代码来实现的话，其实并不容易。线性代数子模块 linalg 提供了一个通用且高效的解决方案，只要把各个方程的常数项写在等号右边，提取出系数数组和常数数组，调用 linalg.solve( ) 函数即可一步求解。

$$\{\begin{array}{r}x-2y+z=0\\ 2y-8z=8\\ -4x+5y+9z=-9\end{array}$$

对于上面的线性方程组，下面的代码演示了求解过程。

>>> A = np.array([[1,-2,1],[0,2,-8],[-4,5,9]]) # 系数数组
>>> b = np.array([0,8,-9]) #常数数组
>>> sla.solve(A, b) # 调用scipy.linalg的solve()函数，返回x、y、z的方程解
array([29., 16., 3.])
>>> nla.solve(A, b) # 调用numpy.linalg的solve()函数，返回x、y、z的方程解
array([29., 16., 3.])