zoukankan      html  css  js  c++  java
  • 使用numpy教你复习线性代数基础

    楔子

    下面我们来一起复习一下线性代数的基础知识,并同时使用numpy进行演示,所以需要你有一些关于numpy的知识(但不需要太多)。另外在线性代数中,存在行列式和矩阵,它们长得都差不多,都类似于二维表的格式。但是行列式要求其行数和列数必须相等,但是矩阵则没有此要求,而我们在创建在numpy中创建行列式和矩阵的时候均使用numpy.array这个函数,这个函数创建的是数组,我们使用数组来模拟行列式和矩阵。这一点要知晓。

    行列式

    全排列与逆序数

    全排列:

    把n个不同的元素排成一列,叫做这n个元素的全排列(或排列)

    n个不同元素进行排列,各个元素的位置不同,那么排列的结果也不同。我们将其所有可能出现的排列的总个数记作(P_{n}),且(P_{n} = n!)。很好理解,n个元素对应n个位置,第一个元素有n种选择、第二个元素有n-1种选择、...、最后一个元素只有一种选择,因此(P_{n})就是n的阶乘。

    import itertools
    import math
    
    
    # 计算n个元素所有可能出现的排列的总个数,以及n的阶乘
    def calc(lst: list):
        return len(list(itertools.permutations(lst, len(lst)))), math.factorial(len(lst))
    
    
    for _ in range(1, 11):
        print(calc(list(range(_))))
    """
    (1, 1)
    (2, 2)
    (6, 6)
    (24, 24)
    (120, 120)
    (720, 720)
    (5040, 5040)
    (40320, 40320)
    (362880, 362880)
    (3628800, 3628800)
    """
    # 手动将n个元素按照不同的顺序进行全排列,然后计算所有可能出现的排列的总个数
    # 然后再计算n的阶乘,当n取1到10的时候,我们验证它们是一样的。
    # 当然无论n取多大它们都是一样的,只不过元素过多的话,创建其所有可能出现的排列会花不少时间
    

    逆序数:

    在一个排列((i_{1},i_{2},...,i_{t},...,i_{s},...i_{n}))中,如果(i_{t}>i_{s}),则称这两个数构成一个逆序。一个排列中所有的逆序的总数称之为此排列的逆序数

    看定义似乎有点难理解,说人话就是:对于排列中的任意一个数,只要它前面存在比它大的数,那么前面那个大的数和该数就构成了一个逆序。

    比如:[3, 2, 4, 1, 5]这个排列:

    第一个元素3前面没有数,所以没有元素能和它构成逆序,那么3这个元素的逆序数是0

    而元素2的前面是3,3比2大,所以3和2就构成了一个逆序,那么2这个元素的逆序数就是1,因为存在一个逆序嘛。

    元素4的前面是3和2,没有比它大的数,所以也没有元素能和它构成逆序,那么4这个元素的逆序数也是0

    元素1的前面是3、2、4,三个元素都比它大,所以都可以和1构成逆序,那么1的逆序数就是3,因为前面存在3个比它大的元素、也就是有3个逆序

    最后一个元素是5,前面没有比它大的元素,所以5的逆序数是0

    然后将每一个元素的逆序数加起来,就称之为这个排列的逆序数,所以上面那个排列的逆序数就是:0 + 1 + 0 + 3 + 0 = 4

    import numpy as np
    
    
    # 计算排列的逆序数
    def calc(lst):
        arr = np.array(lst)
        return np.sum([np.sum(arr[: idx] > val) for idx, val in enumerate(lst)])
    
    # 我们刚刚举例的那个排列
    print(calc([3, 2, 4, 1, 5]))  # 4
    
    # 顺序排列,显然逆序数是0
    print(calc(list(range(10))))  # 0
    # 倒序排列,n个元素的排列,那么逆序数显然是:1 + 2 + 3 + ... + n - 1
    # 在纸上一写马上就能看出,当然聪明如你肯定在脑海里就能明白
    print(calc(list(range(1, 101)[:: -1])))  # 4950 
    

    另外,当一个排列的逆序数为奇数,那么称这个排列为奇排列;当一个排列的逆序数为偶数,那么称这个排列为偶排列

    n阶行列式的概念

    什么是n阶行列式呢?形如:

    (egin{vmatrix}a_{11} & a_{12} & cdots & a_{1n} \ a_{21} & a_{22} & cdots & a_{2n} \vdots & vdots & ddots & vdots\ a_{n1} & a_{n2} & cdots & a_{nn} end{vmatrix})

    这样的式子,我们称之为行列式,它的行数和列数是相等的。其中:(a_{11}、a_{22}、a_{33}、...、a_{nn})成为行列式的主对角线;(a_{1n}、...、a_{n1})称为行列式的次对角线

    另外这个行列式是可以计算的,任何一个行列式,最终计算的结果都是一个标量。但是怎么计算我们后面会说,也会演示如何使用numpy计算一个行列式的值。

    另外,如果把行列式中的某一个元素所在的行和列都去掉,那么会得到一个n-1阶的行列式,这个n-1阶行列式就叫做这个元素的余子式。

    以图中的(a_{23})为例,将其所在的行和列去掉,那么剩下的行列式(阶数显然会减1)就是(a_{23})的余子式。

    因此对于n阶行列式来说,如果将元素(a_{ij})(i表示第几行、j表示第几列)所在的行和列去掉,那么剩下的n-1阶行列式就称之为(a_{ij})的余子式,记作(M_{ij})

    另外:(A_{ij} = (-1)^{i+j}M_{ij}),记作(A_{ij})为元素(a_{ij})的代数余子式

    以上就是余子式和代数余子式。

    n阶行列式的计算

    下面我们来看看如何计算一个n阶行列式,我们先以2阶行列式为例。

    (egin{vmatrix}a & b \ c & d end{vmatrix} = ad - bc)

    所以:(egin{vmatrix}5 & 3 \ 2 & 6 end{vmatrix} = 5 * 6 - 3 * 2 = 24)

    import numpy as np
    
    determinant = np.array([[5, 3], [2, 6]])
    print(determinant)
    """
    [[5 3]
     [2 6]]
    """
    
    # 计算行列式
    print(np.linalg.det(determinant))  # 23.999999999999993
    

    numpy计算行列式有一套自己快速的方法,计算任意阶的行列式都可以使用这个方法。而numpy内部在计算的时候,可能会出现浮点数上的误差,当然这不是numpy的问题,而是计算机在存储浮点数的时候很容易出现这个结果。比如:在计算3和4的平方和、然后再开方,我们都知道是5,但是对于程序来讲,它得出的结果可能是4.9999999999(当然未必会出现,只是举个例子)。总之numpy计算的结果和24在大小上是没有什么区别的,可以认为是相等的。

    print(np.linalg.det(determinant).round())  # 24.0
    

    或者加上一个.round(),会四舍五入取整。

    三阶行列式的计算:

    (egin{vmatrix}a & b & c \ d & e & f \ x & y & z end{vmatrix} = aez + cdy + xbf - cex - zbd - afy)

    所以:(egin{vmatrix}1 & 2 & 3 \ 4 & 5 & 6 \ 7 & 8 & 9 end{vmatrix} = 1 * 5 * 9 + 3 * 4 * 8 + 7 * 2 * 6 - 3 * 5 * 7 - 9 * 2 * 4 - 1 * 6 * 8 = 45 + 96 + 84 - 105 - 72 - 48 = 0)

    import numpy as np
    
    determinant = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
    print(determinant)
    """
    [[1 2 3]
     [4 5 6]
     [7 8 9]]
    """
    print(np.linalg.det(determinant))  # 0.0
    

    四阶行列式的计算:

    四阶行列式的计算,选取任意的一行或者一列,注意:是任意的一行或者一列。然后让每个元素乘上其对应的代数余子式,然后加起来就得到了该四阶行列式的值。

    四阶行列式有四行、四列,随便选择一行或者一列,这里我们选择第一行。然后让每个元素乘上对应的代数余子式,再加起来就是这个四阶行列式的值。我们看到第二项和第四项前面是符号,这是因为2所在的行是第一行、所在的列是第三列,所以((-1)^{1+2} = -1),同理4在第一行、第四列,因此((-1)^{1+4} = -1)。所以,我们说乘的是代数余子式,而不是余子式。

    至于计算的话,我们就不演算了,同理对于更高阶的行列式也是如此。但是显然这个计算量是非常大的,但是不要紧,我们有numpy,numpy内部的方法可以灰常快速的计算行列式的值。

    import numpy as np
    
    determinant = np.arange(1, 17).reshape((4, 4))
    print(determinant)
    """
    [[ 1  2  3  4]
     [ 5  6  7  8]
     [ 9 10 11 12]
     [13 14 15 16]]
    """
    print(np.linalg.det(determinant))  # 0.0
    

    这里我们只介绍到四阶,对于高阶也可以这么做,但是会很麻烦。所以行列式的展开是有技巧的,我们可以通过对行和列进行一些变换来简化运算过程,后面也会说一些行列式在变换时的特点。

    但只会进行说明,至于一个高阶行列式具体应该怎么转换、比如:如何转换成上三角、下三角,以及高斯消元等等,我们只会提一嘴,但不会具体演示怎么计算的。而是会使用numpy,因为笔者是个程序猿,不是学生了,所以对于一个高阶行列式,不会手动计算,而是会直接将数据交给numpy来计算,这一点还请理解。

    import numpy as np
    
    determinant = np.random.randint(1, 10, (8, 8))
    print(determinant)
    """
    [[5 3 7 7 2 6 6 1]
     [6 2 3 5 9 5 3 7]
     [3 1 5 6 5 7 3 7]
     [9 8 9 2 9 1 5 2]
     [7 9 6 4 3 8 1 9]
     [8 2 7 3 3 1 8 7]
     [1 5 6 5 4 4 4 2]
     [3 4 8 5 5 9 6 1]]
    """
    print(np.linalg.det(determinant).round())  # -2904450.0
    

    我们即便是高阶行列式,也是可以使用这个方法计算的,而且速度非常快。即使是100阶的行列式,也可以瞬间出结果。

    n阶行列式的性质

    行列式都具有哪些性质呢?我们来看一下

    1.行列式和它的转置行列式是相等的

    import numpy as np
    
    determinant = np.random.randint(1, 10, (5, 5))
    print(determinant)
    """
    [[4 1 1 7 1]
     [1 9 1 7 9]
     [8 1 1 4 1]
     [8 1 4 8 3]
     [8 4 4 6 2]]
    """
    print(determinant.T)
    """
    [[4 1 8 8 8]
     [1 9 1 1 4]
     [1 1 1 4 4]
     [7 7 4 8 6]
     [1 9 1 3 2]]
    """
    print(np.linalg.det(determinant).round(), np.linalg.det(determinant.T).round())  # -3510.0 -3510.0
    

    2.将行列式的某两行或者某两列交换位置,那么行列式会变号

    import numpy as np
    
    determinant = np.random.randint(1, 10, (5, 5))
    print(determinant)
    """
    [[2 1 9 9 9]
     [2 2 2 1 5]
     [9 1 8 4 9]
     [9 7 3 2 7]
     [3 7 1 2 8]]
    """
    determinant_1 = determinant.copy()
    # 将第1行和第4行交换位置
    determinant_1[[0, 3]] = determinant_1[[3, 0]]
    print(determinant_1)
    """
    [[9 7 3 2 7]
     [2 2 2 1 5]
     [9 1 8 4 9]
     [2 1 9 9 9]
     [3 7 1 2 8]]
    """
    print(np.linalg.det(determinant).round(), np.linalg.det(determinant_1).round())  # -641.0 641.0
    

    3.行列式如果有两行或者两列完全相同,那么此行列式的值为0

    import numpy as np
    
    determinant = np.random.randint(1, 10, (5, 5))
    print(determinant)
    """
    [[7 8 2 5 2]
     [3 4 5 2 5]
     [3 5 1 5 8]
     [6 6 9 8 4]
     [1 5 5 3 1]]
    """
    # 让第3行和第1行相同
    determinant[2] = determinant[0]
    print(determinant)
    """
    [[7 8 2 5 2]
     [3 4 5 2 5]
     [7 8 2 5 2]
     [6 6 9 8 4]
     [1 5 5 3 1]]
    """
    print(np.linalg.det(determinant).round())  # 0.0
    

    4.行列式整体乘上一个常数k,等价于行列式的某一行或者某一列乘上一个常数k

    5.由第4个性质,我们也可以推出,如果某一行或者某一列有公约数,那么这个公约数可以提到行列式的外面。

    6.由第3、第5两个性质,我们也可以推出,如果有两行或者两列成比例,那么这个行列式的值为0

    7.若行列式的某一行或者某一列都是两个元素之和,则该行列式等于两个行列式之和

    我们来解释一下第七个性质:

    import numpy as np
    
    determinant = np.array([
        [1, 2, 3],
        [2, 3, 4],
        [11 + 12, 5 + 8, 17 + 1]
    ])
    
    determinant1 = np.array([
        [1, 2, 3],
        [2, 3, 4],
        [11, 5, 17]
    ])
    
    determinant2 = np.array([
        [1, 2, 3],
        [2, 3, 4],
        [12, 8, 1]
    ])
    
    # 我们将determinant拆分成了determinant1和determinant2
    print(np.linalg.det(determinant).round())  # -15.0
    print(np.linalg.det(determinant1).round(), np.linalg.det(determinant2).round())  # -18.0 3.0
    

    8.将行列式的某一行或者某一列乘上一个常数k,然后加到另外的一行或者一列,行列式的值不变。

    为什么会有这个性质,可以结合前面几条思考一下。

    上三角和下三角行列式

    主对角线下方元素全为0的行列式称为上三角行列式,上方元素全为0的行列式称为下三角行列式。无论是上三角行列式还是下三角行列式,它们的值都等于主对角线上所有元素的乘积。

    import numpy as np
    
    determinant = np.array([
        [1, 0, 0, 0],
        [2, 3, 0, 0],
        [4, 5, 6, 0],
        [4, 5, 6, 8],
    ])
    
    # 我们将determinant拆分成了determinant1和determinant2
    print(np.linalg.det(determinant).round(), 1 * 3 * 6 * 8)  # 144.0 144
    print(determinant.T)
    """
    [[1 2 4 4]
     [0 3 5 5]
     [0 0 6 6]
     [0 0 0 8]]
    """
    print(np.linalg.det(determinant.T).round(), 1 * 3 * 6 * 8)  # 144.0 144
    

    克莱默法则

    如果线性方程组

    (egin{cases}a_{11}x_{1} + a_{12}x_{2} + cdots + a_{1n}x_{n} = b_{1}, \a_{21}x_{1} + a_{22}x_{2} + cdots + a_{2n}x_{n} = b_{2}, \ ... ... ... ... ...\ a_{n1}x_{1} + a_{n2}x_{2} + cdots + a_{nn}x_{n} = b_{n},end{cases})

    的系数行列式D不等于0,那么它有唯一解,(x_{j} = frac{D_{j}}{D}, j=1,2,3,...,n)。其中(D_{j})是把系数行列式D中第j列换成(b_{1},b_{2},b_{3},...,b{n})所得到的行列式。

    如果线性方程组

    (egin{cases}a_{11}x_{1} + a_{12}x_{2} + cdots + a_{1n}x_{n} = 0, \a_{21}x_{1} + a_{22}x_{2} + cdots + a_{2n}x_{n} = 0, \ ... ... ... ... ...\ a_{n1}x_{1} + a_{n2}x_{2} + cdots + a_{nn}x_{n} = 0,end{cases})

    的系数行列式D不等于0,那么它没有非零解。首先根据上面我们知道有唯一解,但如果常数项(b_{1},b_{2},...,b_{n})全是0,那么这个唯一解一定是零解,也就是(x_{1},x_{2},...,x_{n})只可能全部是0。但如果该线性方程组真的存在非零解,那么系数行列式的值必为0,如果不为0,那么不会存在非零解。

    矩阵

    矩阵,我们在numpy中依旧使用array的方式来创建

    矩阵的概念

    由m × n个数(a_{ij}(i=1,2,3,...,m;j=1,2,3,...,n))排列而成的m行n列的数表

    (egin{bmatrix}a_{11}&a_{12}&cdots&a_{1n}\a_{21}&a_{22} &cdots&a_{2n}\vdots&vdots&ddots&vdots\a_{m1}&a_{m2}&cdots&a_{mn} end{bmatrix})

    成为m行n列矩阵,简称m × n矩阵。这m×n个数称之为该矩阵的元素,简称为"元"

    同型矩阵和矩阵相等:

    • 1.两个矩阵的行数和列数都相等,那么称这两个矩阵为同型矩阵
    • 2.两个矩阵为同型矩阵,并且对应位置的元素都相等,那么称这两个矩阵相等。

    因此两个矩阵如果想要相等,那么它的必要条件是这两个矩阵首先要是同型矩阵,否则元素不会一一对应。

    import numpy as np
    
    arr1 = np.array([
        [1, 2, 3],
        [2, 3, 4]
    ])
    
    arr2 = np.array([
        [1, 2, 3],
        [2, 3, 4]
    ])
    
    # arr1和arr2就是同型矩阵
    print(arr1 == arr2)
    """
    [[ True  True  True]
     [ True  True  True]]
    """
    # 但是numpy会返回每一个元素比较的结果
    # 可以使用np.all
    print(np.all(arr1 == arr2))  # True
    

    几种特殊的矩阵

    1.通过矩阵的概念我们得知,矩阵和行列式不同,行列式要求行数和列数必须一致,但是矩阵则没有此要求。但是,如果矩阵的行数和列数也相等的话,那么这个矩阵称之为方阵。行数和列数相等都为n,那么就称这个矩阵为n阶方阵

    2.如果这个矩阵只有一行,称为行矩阵(或行向量);该矩阵只有一列,那么称为列矩阵(或列向量)。所以我们看到向量实际上是矩阵只有一行或者一列的特殊情况。

    但是在numpy中,即使是行矩阵和列矩阵,我们依旧使用二维数组表示

    import numpy as np
    
    arr1 = np.array([[1, 2, 3]])
    arr2 = np.array([[1], [2], [3]])
    print(arr1)
    """
    [[1 2 3]]
    """
    print(arr2)
    """
    [[1]
     [2]
     [3]]
    """
    

    3.主对角线上的元素不全为零,而其它元素全为零的矩阵,称为对角矩阵。记作:(A = diag(λ_{1},λ_{2},...,λ{n}))

    4.元素全部为0的矩阵称之为零矩阵。注意:阶数不同的零矩阵是不相等的,因为它们不是同型矩阵

    5.主对角线上的元素全部为1的对角矩阵,称之为单位矩阵。所以首先是对角矩阵,然后主对角线上的元素全部是1。单位矩阵,一般用E表示

    import numpy as np
    
    
    # numpy生成单位矩阵
    # 表示生成5阶单位矩阵
    arr = np.eye(5)
    print(arr)
    """
    [[1. 0. 0. 0. 0.]
     [0. 1. 0. 0. 0.]
     [0. 0. 1. 0. 0.]
     [0. 0. 0. 1. 0.]
     [0. 0. 0. 0. 1.]]
    """
    
    # 另外
    arr = np.eye(5, M=3)
    # 我们看到如果加上一个参数M=3,那么等于在原来的基础上只选三列
    # M默认为None,如果给定一个具体的值,那么会选择指定个数的列。如果超过了总列数,会用0填充
    print(arr)
    """
    [[1. 0. 0.]
     [0. 1. 0.]
     [0. 0. 1.]
     [0. 0. 0.]
     [0. 0. 0.]]
    """
    
    arr = np.eye(5, k=1)
    # 单位矩阵是主对角线上全部为1
    # 如果指定k=1,那么会在主对角线往上平移一个单位,将对应的元素变成1
    print(arr)
    """
    [[0. 1. 0. 0. 0.]
     [0. 0. 1. 0. 0.]
     [0. 0. 0. 1. 0.]
     [0. 0. 0. 0. 1.]
     [0. 0. 0. 0. 0.]]
    """
    # 同理k为负数,则是往下,这里k=-2
    arr = np.eye(5, k=-2)
    print(arr)
    """
    [[0. 0. 0. 0. 0.]
     [0. 0. 0. 0. 0.]
     [1. 0. 0. 0. 0.]
     [0. 1. 0. 0. 0.]
     [0. 0. 1. 0. 0.]]
    """
    

    另外,一个矩阵的主对角线上所有元素的和称之为矩阵的"迹"

    import numpy as np
    
    arr = np.random.randint(1, 10, (5, 5))
    print(arr)
    """
    [[5 2 9 3 7]
     [1 4 9 6 2]
     [6 8 1 9 9]
     [4 3 4 6 5]
     [9 2 5 7 3]]
    """
    # 计算矩阵的迹
    print(np.trace(arr))  # 19
    

    6.对于一个n阶方阵A,如果满足(a_{ij} = a_{ji}),那么称A为对称矩阵。也就是说,元素关于主对角线是对称的。

    该矩阵:(A = egin{bmatrix}12&6&1\6&8&0\1&0&6end{bmatrix})就是一个对称矩阵,因为元素关于主对角线是对称的。

    7.对于一个矩阵:

    (A = egin{bmatrix}a_{11}&a_{12}&cdots&a_{1n}\a_{21}&a_{22} &cdots&a_{2n}\vdots&vdots&ddots&vdots\a_{m1}&a_{m2}&cdots&a_{mn} end{bmatrix})

    其行列式|A|的每个元素的代数余子式(A_{ij})所构成矩阵

    (A^{*} = egin{bmatrix}A_{11}&A_{12}&cdots&A_{1n}\A_{21}&A_{22} &cdots&A_{2n}\vdots&vdots&ddots&vdots\A_{m1}&A_{m2}&cdots&A_{mn} end{bmatrix})

    称之为A的伴随矩阵

    对于矩阵A和伴随矩阵(A^{*}),有:(AA^{*} = A^{*}A = |A|E)

    矩阵的运算

    矩阵可以做加法,两个同型矩阵相加,等于将矩阵中的对应元素相加。

    import numpy as np
    
    A = np.array([[1, 2, 3], [2, 3, 4]])
    B = np.array([[11, 22, 33], [22, 33, 44]])
    
    print(A)
    """
    [[1 2 3]
     [2 3 4]]
    """
    print(B)
    """
    [[11 22 33]
     [22 33 44]]
    """
    print(A + B)
    """
    [[12 24 36]
     [24 36 48]]
    """
    

    但是在numpy中,即使形状不同也是可以相加的,会进行广播。矩阵相加有如下性质:

    • A + B = B + A
    • (A + B) + C = A + (B + C)
    • 矩阵A里面的所有元素都乘以-1,得到的矩阵成为A的负矩阵,记作-A
    • A + (-A) = 零矩阵,A - B = A + (-B)

    矩阵可以和一个数相乘,等于矩阵内每一个元素都和这个数相乘。

    所以这一点和行列式不同,行列式乘上一个数,等于行列式里面的某一行或者某一列和这个数相乘;但是矩阵和一个数相乘,等于矩阵内每一个元素都和这个数相乘。

    import numpy as np
    
    A = np.array([[1, 2, 3], [2, 3, 4]])
    B = np.array([[11, 22, 33], [22, 33, 44]])
    
    print(A)
    """
    [[1 2 3]
     [2 3 4]]
    """
    print(A * 3)
    """
    [[ 3  6  9]
     [ 6  9 12]]
    """
    

    矩阵和一个数相乘有如下特点:

    • (λμ)A = λ(μA)
    • (λ + μ)A = λA + μA
    • λ(A + B) = λA + λB

    矩阵相加和数乘结合起来,成为矩阵的线性运算。另外在numpy中,支持的运算操作远比目前介绍的丰富。

    矩阵可以和矩阵相乘

    关于矩阵和矩阵相乘,比如:A和B相乘,要求A矩阵的列数必须等于B矩阵的行数。然后A矩阵的每一行乘以B矩阵的每一列,最终得到的矩阵的行数等于A矩阵的行数、列数等于B矩阵的列数

    import numpy as np
    
    A = np.array([[1, 2], [2, 3], [3, 4]])
    B = np.array([[1, 2, 3], [2, 3, 4]])
    print(A)
    """
    [[1 2]
     [2 3]
     [3 4]]
    """
    print(B)
    """
    [[1 2 3]
     [2 3 4]]
    """
    # 注意:这里我们要说明一下,我们目前创建的是数组,不是矩阵。只是用数组来模拟矩阵
    # 所以不能A * B,而是要通过A @ B,来表示矩阵之间的乘法运算
    print(A @ B)
    """
    [[ 5  8 11]
     [ 8 13 18]
     [11 18 25]]
    """
    # 因为*表示乘法运算,如果使用*,那么要先将数组转化成矩阵。
    # 但是np.matrix这个类要被移除了,不推荐使用了
    print(np.matrix(A) * np.matrix(B))
    """
    [[ 5  8 11]
     [ 8 13 18]
     [11 18 25]]
    """
    # 因此虽然我们创建的是数组,但是可以使用@来进行矩阵之间乘法
    # 既然A @ B表示矩阵之间的乘法,那A * B呢?
    # A * B则是相当于矩阵的点乘,也就是A和B对应位置的元素直接相乘
    A = np.array([[1, 2], [2, 3], [3, 4]])
    B = np.array([[1, 2], [2, 3], [3, 4]])
    print(A * B)
    """
    [[ 1  4]
     [ 4  9]
     [ 9 16]]
    """
    # 我们看到就是A和B对应位置的元素进行相乘。
    # 所以矩阵的乘法和矩阵的点乘之间是不同的,矩阵的乘法是A矩阵的每一行和B矩阵的每一列进行内积
    # 但是矩阵的点乘,则是对应位置的元素进行相乘。只不过我们目前创建的是数组,不是矩阵罢了
    # 所以需要使用np.matrix转成矩阵,再相乘。或者数组之间使用@操作符,表示直接按照矩阵相乘的方式进行运算。
    # 至于数组本身之间的乘法,则是类似矩阵的点乘,也就是对应位置的元素进行相乘
    

    矩阵相乘具有如下性质:

    • 矩阵乘法一般不满足交换律,即AB≠BA。其实一般来说,AB可以相乘,往往BA并不能相乘
    • 若矩阵满足AB=零矩阵,并不能说明A或者B中间有一个是零矩阵
    import numpy as np
    
    A = np.array([[1, -1], [1, -1]])
    B = np.array([[1, -1], [1, -1]])
    print(A @ B)
    """
    [[0 0]
     [0 0]]
    """
    
    • (AB)C = A(BC)
    import numpy as np
    
    A = np.array([[2, 3], [1, -7]])
    B = np.array([[-5, 1], [2, -5]])
    C = np.array([[4, 2], [2, 6]])
    print((A @ B) @ C)
    """
    [[-42 -86]
     [ -4 178]]
    """
    print(A @ (B @ C))
    """
    [[-42 -86]
     [ -4 178]]
    """
    
    • A(B + C) = AB + AC
    import numpy as np
    
    A = np.array([[2, 3], [1, -7]])
    B = np.array([[-5, 1], [2, -5]])
    C = np.array([[4, 2], [2, 6]])
    print(A @ (B + C))
    """
    [[ 10   9]
     [-29  -4]]
    """
    print(A @ B + A @ C)
    """
    [[ 10   9]
     [-29  -4]]
    """
    
    • λ(AB) = (λA)B = A(λB)
    import numpy as np
    
    A = np.array([[2, 3], [1, -7]])
    B = np.array([[-5, 1], [2, -5]])
    
    print(3 * (A @ B))
    """
    [[-12 -39]
     [-57 108]]
    """
    print((3 * A) @ B)
    """
    [[-12 -39]
     [-57 108]]
    """
    print(A @ (3 * B))
    """
    [[-12 -39]
     [-57 108]]
    """
    
    • AE = EA = A, 单位矩阵E在矩阵乘法中类似于数字运算中的1
    import numpy as np
    
    A = np.array([[2, 3], [1, -7]])
    print(A)
    """
    [[ 2  3]
     [ 1 -7]]
    """
    print(A @ np.eye(2))
    """
    [[ 2.  3.]
     [ 1. -7.]]
    """
    
    • 若A是n阶方阵,则A的m次方乘上A的k次方,等于A的m+k次方;A先m次方、然后整体在k次方,等于A的m乘k次方。注意:但是一般来说,AB的k次方,不等于A次方乘上B的k次方
    import numpy as np
    
    A = np.random.randint(1, 10, (4, 4))
    # 我们说数组之间的@,等价于矩阵的*
    # 但是数组来说,如果使用@实现矩阵的幂次方运算比较费劲,其实也不算费劲
    # 不过最方便的办法还是将这里的数组转为矩阵吧,这样就可以使用**操作符了
    A = np.matrix(A)
    print(np.all(A ** 5 * A ** 8 == A ** (5 + 8)))  # True
    
    print(np.all(
        (A ** 5) ** 8 == A ** (5 * 8)
    ))  # True 
    

    矩阵的转置

    把矩阵的A的行换成同序数的列所得到的矩阵叫做A的转置矩阵,记作:(A^{T})

    import numpy as np
    
    A = np.random.randint(1, 10, (2, 3))
    print(A)
    """
    [[6 1 9]
     [9 4 8]]
    """
    # 虽然我们创建的是数组,但是大部分情况下可以直接当成矩阵来用
    # 直接调用A.T即可得到转置矩阵
    print(A.T)
    """
    [[6 9]
     [1 4]
     [9 8]]
    """
    

    注意:如果A为对称矩阵,那么有(A = A^{T})

    转置矩阵有如下性质:

    • ((A^{T})^{T} = A)
    import numpy as np
    
    A = np.random.randint(1, 10, (2, 3))
    print(A)
    """
    [[7 4 3]
     [9 9 9]]
    """
    print(A.T.T)
    """
    [[7 4 3]
     [9 9 9]]
    """
    
    • ((A + B)^{T} = A^{T} + B^{T})
    import numpy as np
    
    A = np.random.randint(1, 10, (2, 3))
    B = np.random.randint(1, 10, (2, 3))
    print((A + B).T)
    """
    [[ 8  7]
     [15 12]
     [13 11]]
    """
    print(A.T + B.T)
    """
    [[ 8  7]
     [15 12]
     [13 11]]
    """
    
    • ((λA)^{T} = λA^{T})
    import numpy as np
    
    A = np.random.randint(1, 10, (2, 3))
    print((3 * A).T)
    """
    [[18  3]
     [24 15]
     [ 9 21]]
    """
    print(3 * A.T)
    """
    [[18  3]
     [24 15]
     [ 9 21]]
    """
    
    • ((AB)^{T} = B^{T}A^{T})
    import numpy as np
    
    A = np.random.randint(1, 10, (2, 3))
    B = np.random.randint(1, 10, (2, 3))
    print((A * B).T)
    """
    [[18  2]
     [49  6]
     [ 2 18]]
    """
    print(B.T * A.T)
    """
    [[18  2]
     [49  6]
     [ 2 18]]
    """
    

    方阵的行列式

    由n阶方阵A的所有元素所构成的行列式,叫做方阵A的行列式,记作|A|或者det A。如果该方阵的行列式值为0,那么称这个矩阵为奇异矩阵;值不为0,那么称这个矩阵为非奇异矩阵。

    • (|A^{T}| = |A|)
    import numpy as np
    
    A = np.random.randint(1, 10, (3, 3))
    print(np.linalg.det(A).round(), np.linalg.det(A.T).round())  # 248.0 248.0
    
    • (|λA| = λ^{n}|A|),这里是λ的n次方,很好理解,我们说矩阵和数相乘,那么这个数会作用于所有的元素,因此相当于每一行都提出来一个λ。总共有n行,那么对于行列式来讲就是λ的n次方。
    import numpy as np
    
    A = np.random.randint(1, 10, (3, 3))
    print(np.linalg.det(A * 5).round(), np.linalg.det(A).round() * 5 ** 3)  # 12125.0 12125.0
    
    • (|AB| = |A||B|)
    import numpy as np
    
    A = np.random.randint(1, 10, (3, 3))
    B = np.random.randint(1, 10, (3, 3))
    print(np.linalg.det(A @ B).round(), (np.linalg.det(A) * np.linalg.det(B)).round())  # 3737.0 3737.0
    
    • (|AB| = |BA|)
    import numpy as np
    
    A = np.random.randint(1, 10, (3, 3))
    B = np.random.randint(1, 10, (3, 3))
    print(np.linalg.det(A @ B).round(), np.linalg.det(B @ A).round())  # 2622.0 2622.0
    

    逆矩阵

    对于n阶方阵A,如果有一个n阶方阵B,是的AB=BA=E,那么说方阵A是可逆的,并把方阵B叫做方阵A的逆矩阵。另外,A的逆矩阵也叫作(A^{-1})

    import numpy as np
    
    A = np.random.randint(1, 10, (3, 3))
    # 求逆矩阵B
    B = np.linalg.inv(A)
    print(B)
    """
    [[ 0.25       -0.22916667  0.14583333]
     [-0.75        0.60416667 -0.02083333]
     [ 1.         -0.58333333 -0.08333333]]
    """
    print((A @ B).round())
    """
    [[ 1.  0. -0.]
     [-0.  1. -0.]
     [ 0.  0.  1.]]
    """
    

    我们说一个矩阵如果有逆矩阵,那么这个矩阵首先要是一个方阵。但是在numpy中,即使不是方阵, 也可以求出逆矩阵。

    import numpy as np
    
    A = np.random.randint(1, 10, (3, 5))
    # 求逆矩阵使用inv,如果不是方阵,那么使用pinv
    B = np.linalg.pinv(A)
    print(B)
    """
    [[-0.07192369  0.07431035  0.04447579]
     [ 0.01839564 -0.00424631  0.02602996]
     [ 0.05700139 -0.08014275  0.07874739]
     [-0.01293216 -0.00438183  0.037499  ]
     [ 0.05127439  0.07520378 -0.105209  ]]
    """
    print((A @ B).round())
    """
    [[ 1.  0. -0.]
     [-0.  1. -0.]
     [ 0.  0.  1.]]
    """
    

    逆矩阵有如下性质:

    • 方阵A可逆的充要条件是|A|≠0,且(A^{-1} = frac{1}{|A|}A^{*}),其中(A^{*})为A的伴随矩阵
    import numpy as np
    
    A = np.array(
        [
            [2, 2],
            [3, 3]
        ]
    )
    # 此矩阵对应的行列式值为0
    print(np.linalg.det(A).round())  # 0.0
    # 我们来求A的逆矩阵,我们说一个矩阵如果有逆矩阵,那么行列式的值必不为0
    try:
        print(np.linalg.inv(A))
    except Exception as e:
        print(e)  # Singular matrix
    # 报错了,提示这个矩阵是奇异矩阵
    
    • 若方阵A可逆,则(A^{-1})可逆,且((A^{-1})^{-1} = A)
    import numpy as np
    
    A = np.array(
        [
            [5, 12],
            [30, 3]
        ]
    )
    print(np.linalg.inv(np.linalg.inv(A)))
    """
    [[ 5. 12.]
     [30.  3.]]
    """
    
    • 若方阵A可逆,且λ≠0,则λA可逆,且((λA)^{-1} = frac{1}{λ}A^{-1})
    import numpy as np
    
    A = np.array(
        [
            [5, 12],
            [30, 3]
        ]
    )
    print(np.linalg.inv(3 * A))
    """
    [[-0.00289855  0.0115942 ]
     [ 0.02898551 -0.00483092]]
    """
    print(np.linalg.inv(A) * 1 / 3)
    """
    [[-0.00289855  0.0115942 ]
     [ 0.02898551 -0.00483092]]
    """
    
    • 若方阵A,B为同型方阵且均可逆,则AB可逆,且((AB)^{-1} = B^{-1}A^{-1})
    import numpy as np
    
    A = np.random.randint(1, 10, (3, 3))
    B = np.random.randint(1, 10, (3, 3))
    print(np.linalg.inv(A @ B))
    """
    [[ 0.01923077 -0.05        0.04807692]
     [-0.00925926  0.08703704 -0.10648148]
     [-0.01210826 -0.00925926  0.05306268]]
    """
    print(np.linalg.inv(B) @ np.linalg.inv(A))
    """
    [[ 0.01923077 -0.05        0.04807692]
     [-0.00925926  0.08703704 -0.10648148]
     [-0.01210826 -0.00925926  0.05306268]]
    """
    
    • 若方阵A可逆,则(A^{T})可逆,且((A^{T})^{-1} = (A^{-1})^{T})

      另外,当|A|≠0、k为整数时,定义(A^{0} = E)(A^{-k} = (A^{-1})^{k});当|A|≠0,λ、μ为整数时,(A^{λ}A^{μ} = A^{λ + μ})((A^{λ})^{μ} = A^{λμ}),当然这个我们上面其实说过了

    import numpy as np
    
    A = np.random.randint(1, 10, (3, 3))
    print(np.linalg.inv(A.T))
    """
    [[ 0.16666667  0.16666667 -0.16666667]
     [-0.04545455 -0.40909091  0.31818182]
     [-0.15151515  0.3030303   0.06060606]]
    """
    print(np.linalg.inv(A).T)
    """
    [[ 0.16666667  0.16666667 -0.16666667]
     [-0.04545455 -0.40909091  0.31818182]
     [-0.15151515  0.3030303   0.06060606]]
    """
    
    • 若A可逆,则有(|A^{-1}| = |A|^{-1})
    import numpy as np
    
    A = np.random.randint(1, 10, (3, 3))
    print(np.linalg.det(np.linalg.inv(A.T)))
    """
    -0.010869565217391302
    """
    print(np.linalg.det(A) ** -1)
    """
    -0.010869565217391302
    """
    

    矩阵的秩

    在m×n矩阵A中任取k行k列(k≤m,k≤n,按照从上到下、从左到右的顺序),所得到的k阶行列式,称之为A的k阶子式。

    而矩阵A的秩就是A中不为0的子式的最高阶数,记作R(A)

    import numpy as np
    
    A = np.random.randint(1, 10, (8, 8))
    print(A)
    """
    [[3 2 9 8 1 4 9 8]
     [9 8 8 5 6 8 9 4]
     [1 9 2 1 2 9 3 7]
     [8 2 5 2 6 1 5 5]
     [1 6 3 1 4 1 4 2]
     [5 5 9 9 2 5 8 5]
     [1 7 9 1 3 5 5 9]
     [1 5 6 2 6 1 1 9]]
    """
    # matrix_rank表示求矩阵的秩
    print(np.linalg.matrix_rank(A))  # 8
    # 将最后一行全部改成0
    A[-1] = 0
    print(np.linalg.matrix_rank(A))  # 7
    
    """
    首先该矩阵A是8阶方阵,那么显然它的子式最高阶就是8,也就是矩阵A本身对应的行列式。而该行列式的值是不为0的,所以它的秩是8
    然后我们将最后一行变成了0,显然A本身对应的行列式的值变成了0,因此不存在不为0的8阶子式。而且8阶子式只有一个,就是A本身对应的行列式
    既然8阶没有,那么选取7阶,显然7阶是成立的。另外如果某个k阶子式有多个,只要找到一个不为0的子式即可。
    """
    # 其实如果你仔细学过线性代数的话,那么你应该知道一种很方便的求秩的方法
    # 那就是将该矩阵通过线性变换转化成阶梯型,然后非零行的个数就是该矩阵的秩
    # 不过至于怎么变换我们这里就不说了,因为计算的时候数据量一般会非常大,肯定不会手算。
    # 我们只需要知道它的含义就可以了,至于计算的工作交给计算机
    

    n维向量

    向量以及向量组的概念

    n个有次序的数(a_{1}、a_{2}、a_{3}、...、a_{n})所组成的数组称为n维向量(但是说实话,这样的一个数组在numpy中实际上是一维的,里面有n个元素,但在数学上定义是n维的,所以要注意),这n个数称为该向量的n个分量,第i个数(a_{i})称为第i个分量。

    若干个同维数的行向量或者列向量所组成的集合叫做向量组

    向量组的线性相关性

    向量组表示向量

    假设存在一个向量β和一个向量组(A=(α_{1},α_{2},...,α_{n})),那么我们可以使用向量组来表示该向量,如:(β=λ_{1}α_{1}λ_{2}α_{2}+...+λ_{n}α_{n})。但是向量β能被向量A表示的充要条件是矩阵(A=(α_{1},α_{2},...,α_{n}))的秩等于矩阵(B=(α_{1},α_{2},...,α_{n},b))的秩

    线性相关和线性无关

    对于一个向量组(A=(α_{1},α_{2},...,α_{n})),如果存在一组不全为0的数(k_{1},k_{2},...,k_{n}),使得(k_{1}α_{1}+k_{2}α_{2}+...+k_{n}α_{n}=vec{0}),则称向量组A是线性相关的,否则称它为线性无关的。

    所以由以上可得出:

    • 任意一个向量组不是线性相关就是线性无关
    • 包含零向量的向量组一定线性相关
    • 单个非零向量组成的向量组一定是线性无关的
    • 两个向量组成的向量组线性相关的充要条件是这两个向量对应的分量成比例

    另外,一个向量组(A=(α_{1},α_{2},...,α_{n}))线性相关的充分必要条件是它所构成的矩阵的秩小于向量组里向量的个数n,即R(A) < m;线性无关的充要条件则是R(A) = m

    向量组的线性相关性有如下性质:

    • 向量组内部分向量相关则整个向量组相关。很好理解,因为是一组不全为0的数,k1,k2...如果部分相关,那么让剩下的全部为0即可。所以部分相关,那么整体也是相关的
    • 线性无关的向量组,将分量延长后也是线性无关的
    • m个n维向量,当维数n小于m时,那么这个m个向量组成的向量组一定是线性相关的
    • 设向量组A: a1,a2,...,an线性无关,向量组B: a1,a2,...,an,b线性相关,则向量b一定可以被向量组A线性表示,并且是唯一的

    最大无关组与向量组的秩

    对于一个向量组A,如果能在A中选出r个向量,(α_{1},α_{2},...,α_{r}),满足向量组(A_{0}: α_{1},α_{2},...,α_{r})线性无关,并且向量组A中任意r+1个向量都线性相关,那么称向量组(A_{0})是向量组A的一个最大线性无关组。并且,最大线性无关组中向量的个数称为向量组A的秩。

    注意:只包含零向量的向量组没有最大线性无关组,规定它的秩为0

    所以向量组和矩阵还是比较类似,都是由多个行向量或者列向量所构成。所以矩阵的秩等于它行向量组的秩,也等于它列向量组的秩。

    另外:如果一个向量组B能由向量组A线性表示,则向量组B的秩不大于向量组A的秩。

    线性方程组

    线程方程组有解的判定条件

    定理一:n元齐次线性方程组(A_{m×n}x=0)有非零解的充要条件是系数矩阵的秩R(A)<n

    定理二:n元非齐次线性方程组(A_{m×n}x=b)有解的充要条件是系数矩阵的秩R(A)等于增广矩阵B=(A, b)的秩R(B)

    线程方程组解的结构

    基础解系的定义

    若向量组(η_{1},η_{2},η_{3},...,η_{t})满足:(η_{1},η_{2},η_{3},...,η_{t})是Ax=0的一组线性无关的解,并且Ax=0的任意解都可由(η_{1},η_{2},η_{3},...,η_{t})线性表示,那么称(η_{1},η_{2},η_{3},...,η_{t})为齐次线性方程组Ax=0的基础解系。

    如果向量组(η_{1},η_{2},η_{3},...,η_{t})齐为次线性方程组Ax=0的基础解系,那么Ax=0的通解可以表示为:(x=k_{1}η_{1}+k_{2}η_{2}+...+k_{n}η_{n}, 其中k_{1},k_{2},...,k_{n}是任意常数)

    但是说实话,我个人觉得一般方程存在基础解系,那么基本上是等式的个数小于未知数的个数。比如:有三个等式,但是存在4个变量。

    线性方程组的解法

    形如(A_{m×n}=0)的齐次线性方程组:系数矩阵化成行最简矩阵,便可写出其通解

    形如(A_{m×n}=b)的非齐次线性方程组:增广矩阵化成阶梯型矩阵,便可判断其是否有解。如果有解, 化成行最简矩阵,便可写出其通解。

    至于怎么化简可以自己去搜索,这里我们来演示一下如何使用numpy来解方程

    import numpy as np
    
    """
    x+2y+z=7
    2x-y+3z=7
    3x+y+2z=18
    求这个方程组的解
    """
    # 首先将系数写下来,排成一个矩阵,也就是得到系数矩阵
    A = np.array([[1, 2, 1],
                  [2, -1, 3],
                  [3, 1, 2]])
    # 将右边的常数写下来,排成一个矩阵
    b = np.array([7, 7, 18])
    # 求解,将参数传进去
    x = np.linalg.solve(A, b)
    print(x)  # [ 7.  1. -2.]
    
    # 验证,将系数矩阵和解相乘,会得到b
    print(A @ x)  # [ 7.  7. 18.]
    

    矩阵的特征值和特征向量

    定义

    假设A是n阶方阵,如果数λ和n维非零列向量x是关系式Ax=λx成立,那么这样的数λ称为方阵A的特征值,非零向量称为A的特征向量。

    所以把两边的x约掉,特征方程|A-λE|=0的根就是方阵A的特征值,注意:x是个向量,约掉得到单位矩阵E,而不是1。齐次线性方程(A - λE)=0的非零解就是对应的特征向量。

    numpy其特征值和特征向量

    import numpy as np
    
    A = np.array([[3, -1], [-1, 3]])
    w, v = np.linalg.eig(A)
    print("特征值:", w)
    """
    特征值: [4. 2.]
    """
    print("特征向量:", v)
    """
    特征向量: [[ 0.70710678  0.70710678]
     [-0.70710678  0.70710678]]
    """
    # 这个特征向量的值不止一种
    # 当特征值为4的时候,特征向量为(x1, x2)且x1-x2=0
    # 当特征值为2的时候,特征向量为(x1, x2)且x1+x2=0
    
  • 相关阅读:
    组合数问题
    [Repost] 悬线法
    图论 List
    杂项 List
    动态规划 List
    Binary Search
    树状数组,Fenwick Tree
    HDU1086判断线段相交
    高效大数模板
    HDUOJ2298三分加二分
  • 原文地址:https://www.cnblogs.com/traditional/p/12658555.html
Copyright © 2011-2022 走看看