为什么numpy在反转奇异矩阵时不会引发错误？

问：

我正在编写一个教程，作为其中的一部分，我想展示在分析具有奇异协方差矩阵的数据时出现的一些问题。

我在numpy中遇到了一些奇怪的行为。

有问题的行为是这样的：

1. 给定一个奇异矩阵 （A），np.linalg.inv（A） 将执行而不会引发异常
2. 该函数返回一个矩阵，该矩阵不是原始矩阵 A*inv（A） ！= I 的倒数

我的问题本质上是：谁能解释为什么允许这种行为，或者这是一个错误，或者我误解了什么？我不是在寻找如何获得奇异矩阵的逆或近似逆的解决方案。

此外，numpy 使用什么方法来测试奇点？它似乎没有单独检查条件编号（见下文）。

这个问题似乎与相同的行为有关（至少上面的第 1 点），但是公认的答案（由于协方差矩阵中的数字非常小，OP 首先测试他们的矩阵是奇异的，因此无法正常工作不适用于我的情况）。

下面是一个完整的示例，演示了该问题。我正在使用 Python 1.24.2 在 numpy 3.10 中工作。任何帮助都非常感谢。

Numpy 无错误地反转奇异矩阵

import numpy as np
from scipy.stats import ortho_group

def generate_random_cov_matrix(sz, rank, eig_val_scale):
    """
    Generate a random covariance matrix with known rank by building it from a specified number of eigenvectors and eigenvalues

    :param sz: the returned covariance matrix with be shaped sz x sz
    :param rank: the desired rank of the returned matrix
    :param eig_val_scale: eigenvalues will be randomly sampled from a chi2 distribution (df=1) and scaled by this parameter
    :return: the covariance matrix

    """
    if rank > sz:
        raise ValueError('rank cannot be greater than size')

    #
    eig_vecs = ortho_group.rvs(dim=sz)
    eig_vecs = eig_vecs[:, 0:rank]
    # randomly sample eigenvalues from chi square (df=1) just so that  smaller eigenvalues are more likely
    eig_vals = np.random.chisquare(1, rank) * eig_val_scale
    # sort in descending order
    eig_vals = np.sort(eig_vals)
    eig_vals = eig_vals[::-1]
    # build matrix from eigenvectors and eigenvalues
    cov = eig_vecs @ np.diag(eig_vals) @ eig_vecs.T
    return cov, eig_vecs, eig_vals


# generate a singular 20 x 20 covariance matrix
n_variables = 20
rank = 10
cov_mat,_,_ = generate_random_cov_matrix(n_variables,rank,100)

# is it singular as expected
assert np.isclose(np.linalg.det(cov_mat),[0.])

# are the elements of the of the covariance matrix so small that their determinant could be >0 but zero due to machine precision
print(np.mean(cov_mat))


# does numpy compute its inverse without an exception
try:
    cov_inv = np.linalg.inv(cov_mat)
    print('No exception was raised')
except:
    print('An exception was raised')


# is the resulting inverse a proper inverse (is A*inv(A) equal to an identity matrix)

try:
    assert np.isclose(cov_mat @ cov_inv,np.identity(n_variables))
except AssertionError:
    print('The computed inverse is not a proper inverse')

Numpy 不会单独检查条件编号

# runs without failure
arr = np.array([[1.0, 1.2], [2.0 + 1e-14, 2.4 + 1e-14]])
rank = np.linalg.matrix_rank(arr)
cond = np.linalg.cond(arr)
np.linalg.inv(arr)
print("rank: ", rank)
print("cond: ", cond)

# fails despite having the same condition number
arr = np.array([[1.0, 1.2], [2.0 + 1e-15, 2.4 + 1e-15]])
rank = np.linalg.matrix_rank(arr)
np.testing.assert_approx_equal(cond, np.linalg.cond(arr))
cond = np.linalg.cond(arr)
print("rank: ", rank)
print("cond: ", cond)
np.linalg.inv(arr)