使用 scipy 计算矩阵排名

Posted 2023-02-23

【中文标题】使用 scipy 计算矩阵排名

【英文标题】：Calculate Matrix Rank using scipy

【发布时间】：2011-01-29 06:15:21

【问题描述】：

我想使用 scipy 计算矩阵的mathematical rank。最明显的函数numpy.rank 计算数组的维数（即标量的维数为 0，向量为 1，矩阵为 2，等等...）。我知道numpy.linalg.lstsq 模块具有此功能，但我想知道这样的基本操作是否内置在某处的矩阵类中。

这是一个明确的例子：

from numpy import matrix, rank
A = matrix([[1,3,7],[2,8,3],[7,8,1]])
print rank(A)

这给出了2 的维度，我正在寻找3 的答案。

【问题讨论】：

我使用 Mathematica 检查了排名 - 确实是 3。您在 Python 中调用的函数要么不正确，要么使用错误。

用法是正确的——这首先让我感到困惑。在帖子中，我解释了 rank 的作用：它计算数组的维数。 “rank-3”数组将是 list-of-lists-of-lists。

请注意，“排名”一词有些含糊。对于张量，等级告诉您索引的数量（例如，标量是等级为 0 的张量，向量等级为 1，矩阵等级为 2）。对于线性代数，还有你上面引用的定义。从文档字符串中，很明显 Numpy 使用了前者。

【参考方案1】：

Numpy 提供numpy.linalg.matrix_rank():

>>> import numpy
>>> numpy.__version__
'1.5.1'
>>> A = numpy.matrix([[1,3,7],[2,8,3],[7,8,1]])
>>> numpy.linalg.matrix_rank(A)
3

【讨论】：

我怎样才能找到 integer 矩阵 modulo n 的等级？在 Mathematica 中有这个函数 MatrixRank[..., Modulus -> n]，但是如何在 Python 中实现这个函数呢？

【参考方案2】：

为需要在实践中完成这项工作的人提供粗略的代码 sn-p。随时改进。

u, s, v = np.linalg.svd(A)
rank = np.sum(s > 1e-10)

【参考方案3】：

如果numpy 不提供排名功能，您为什么不自己编写？

计算秩的一种有效方法是通过奇异值分解 - 矩阵的秩等于非零奇异值的数量。

def rank(A, eps=1e-12):
    u, s, vh = numpy.linalg.svd(A)
    return len([x for x in s if abs(x) > eps])

请注意，eps 取决于您的应用程序 - 大多数人会同意 1e-12 对应于零，但即使 eps=1e-9，您也可能会看到数值不稳定。

使用您的示例，答案是三个。如果将第二行更改为 [2, 6, 14]（与第一行线性相关）答案是二（“零”特征值为 4.9960E-16）

【参考方案4】：

此答案已过时。

答案是否定的——目前在 scipy 中没有专门用于计算数组/矩阵的矩阵秩的函数。添加一个之前已经讨论过，但如果它会发生，我不相信它还没有。

【讨论】：

顺便提一下，mail.scipy.org/pipermail/numpy-discussion/2008-February/… 是关于此问题的简短新闻组讨论的第一篇文章。

现在有numpy.linalg.matrix_rank()。看我的回答。

【参考方案5】：

我不特别了解 Numpy，但这不太可能是矩阵的内置操作；它涉及相当密集的数值计算（以及对浮点舍入误差等的相关关注）和阈值选择，在给定的上下文中可能合适也可能不合适，算法选择对于准确快速地计算很重要。

基本类中内置的东西往往是可以以独特而直接的方式执行的东西，例如最复杂的矩阵乘法。

【讨论】：

这是一个很好的观点，一个数值不稳定的矩阵可能会由于舍入误差而导致排名发生变化。但是，这是一个已知问题，我想知道 scipy/numpy 库是否直接具有功能。如果答案是否定的——那也没关系，我总是可以选择 SVD。

这不仅仅是数值不稳定的。 1.0, 3.0, 1.0/3.0, 1.0 怎么样？分工不能给出准确的答案，所以这应该算作第一名还是第二名？

【参考方案6】：

线性代数函数通常分组在numpy.linalg。 （它们也可以从 scipy.linalg 获得，它具有更多功能。）这允许多态性：函数可以接受 SciPy 处理的任何类型。

所以，是的，numpy.linalg.lstsq 函数可以满足您的要求。为什么不够？

【讨论】：

它可以满足我的要求 - 但它做了很多不必要的事情，并且有大量的行李。同样可以通过 LU 分解然后行排序来完成。问题的意图 - 如果不清楚，是否存在一个函数，其唯一目的是计算排名。 IE。接受一个矩阵，吐出一个整数。

【参考方案7】：

scipy 现在包含一个高效的interpolative method，用于使用随机方法估计矩阵/LinearOperator 的秩，这通常足够准确：

>>> from numpy import matrix
>>> A = matrix([[1,3,7],[2,8,3],[7,8,1]], dtype=float)  # doesn't accept int

>>> import scipy.linalg.interpolative as sli
>>> sli.estimate_rank(A, eps=1e-10)
3