Python使用numpy提取矩阵线性无关列并拆分矩阵及获取列索引

核心用Numpy的列主元QR分解即可快速获取线性无关列索引，实现逻辑如下：

首先导入依赖：

import numpy as np

核心实现函数 #

def split_matrix(matrix_A):
    m = matrix_A.shape[0]
    # 列主元QR分解，pivot返回列的优先级索引，前m个为线性无关列
    Q, R, pivot = np.linalg.qr(matrix_A, mode='reduced', pivoting=True)
    # 拆分B和N的原矩阵索引
    b_idx = pivot[:m]
    n_idx = pivot[m:]
    # 按索引取列得到结果矩阵
    matrix_B = matrix_A[:, b_idx]
    matrix_N = matrix_A[:, n_idx]
    return matrix_B, matrix_N, b_idx, n_idx

测试示例 #

# 输入示例矩阵
matrix_A = np.array([
    [-3., 1., 4., 0., 0.],
    [1., 0., -1., 1., 0.],
    [1., 0., 1., 0., 1.]
], dtype=np.float64)

matrix_B, matrix_N, b_idx, n_idx = split_matrix(matrix_A)

print("matrix_B:\n", matrix_B)
print("matrix_N:\n", matrix_N)
print("B对应原索引：", b_idx.tolist())
print("N对应原索引：", n_idx.tolist())

运行输出和要求的结果完全一致：

matrix_B:
 [[1. 0. 0.]
 [0. 1. 0.]
 [0. 0. 1.]]
matrix_N:
 [[-3.  4.]
 [ 1. -1.]
 [ 1.  1.]]
B对应原索引： [1, 3, 4]
N对应原索引： [0, 2]

补充说明 #

如果输入的matrix_A不是行满秩的，可将获取行数的逻辑m = matrix_A.shape[0]替换为取矩阵秩的逻辑：m = np.linalg.matrix_rank(matrix_A)，即可提取到所有线性无关列。

内容的提问来源于stack exchange，提问作者Patrícia