反转协方差矩阵的马氏距离

Posted 2023-02-19

【中文标题】反转协方差矩阵的马氏距离

【英文标题】：Mahalanobis distance inverting the covariance matrix

【发布时间】：2012-07-31 21:08:53

【问题描述】：

我正在编写一个函数来获取两个向量之间的马氏距离。我知道这是使用方程 a'*C^-1*b 实现的，其中 a 和 b 是向量，C 是协方差矩阵。我的问题是，有没有一种有效的方法可以在不使用 Gauss-Jordan 消除的情况下找到矩阵的逆矩阵，或者没有办法解决这个问题？我正在寻找一种自己做的方法，而不是使用任何预定义的函数。

我知道 C 是 Hermitian 正定矩阵，那么有什么方法可以在算法上利用这个事实吗？还是有一些聪明的方法可以计算马氏距离而不计算协方差的倒数？任何帮助将不胜感激。

***编辑：上面的马氏距离方程不正确。它应该是 x'*C^-1*x 其中 x = (b-a)，b 和 a 是我们试图找到其距离的两个向量（感谢 LRPurser）。因此，所选答案中提出的解决方案如下：

d=x'*b，其中 b = C^-1*x C*b = x，因此使用 LU 分解或 LDL' 分解求解 b。

【问题讨论】：

Gauss-Jordan 消元 是 求逆矩阵（或至少求解线性系统）的有效方法。至少其中之一。

您可以使用Cholesky decomposition 求解具有正定系统矩阵的线性方程组。

【参考方案1】：

您可以（并且应该！）使用 LU decomposition 而不是显式反转矩阵：使用分解求解 C x = b 比计算 C^-1 并乘以向量 b 具有更好的数值属性。

由于您的矩阵是对称的，因此 LU 分解实际上等效于 LDL* decomposition，这是您实际应该使用的。由于您的矩阵也是正定矩阵，因此您应该能够在不旋转的情况下执行此分解。

---

编辑：请注意，对于这个应用程序，您不需要解决完整的C x = b 问题。

相反，给定 C = L D L* 和差分向量 v = a-b，求解 L* y = v 为 y（这是完整 LU 求解器的一半工作量）。

然后，y^t D^-1 y = v^t C^-1 v 可以在线性时间内计算出来。

【讨论】：

所以解系统C*x=b，然后乘a'*b得到范数，就完美了！正是我想要的。非常感谢！

只是出于好奇，你知道 QR 分解吗？我的理解是它在数学上比 LU 分解更稳定，这也是 MATLAB 在求解线性系统时使用的。您知道在这种情况下这与 LDL* 分解相比如何吗？

我不知道 QR 可能会给你带来什么，抱歉。如果您认为您可能在使用未旋转的 LDL* 时存在稳定性问题，那么也请尝试使用 QR 看看是否有帮助...

没问题。我通常使用 QR 分解而不是 LU 分解，但在这种情况下，我怀疑 LDL* 分解实际上可能更好。再次感谢！

【参考方案2】：

第一马氏距离 (MD) 是关于两个向量测量中不确定性的标准化距离。当C=Indentity matrix 时，MD 减少到欧几里得距离，因此乘积减少到向量范数。此外，对于所有非零向量，MD 始终为正定或大于零。 By your formulation with the appropriate choice of the vectors a and b, a*C^-1*b can be less than zero.希望您正在寻找的向量的差异是x=(a-b)，这使得计算x^t*C^-1*x，其中x^t 是向量x 的转置。另请注意MD=sqrt(x^t*C^-1*x) 由于您的矩阵是对称且正定的，因此您可以利用 Cholesky 分解(MatLab-chol)，它使用一半的运算作为 LU，并且在数值上更稳定。 chol(C)=L 其中C=L*L^t 其中L 是下三角矩阵，L^t 是L 的转置，使其成为上三角矩阵。你的算法应该是这样的

（Matlab）

  x=a-b;
  L=chol(C); 
  z=L\x; 
  MD=z'*z; 
  MD=sqrt(MD);

【讨论】：

根据wiki.scipy.org/NumPy_for_Matlab_Users，chol(a) 返回一个上三角矩阵。

MD=z'*z; MD=sqrt(MD); 只是 z 的标准吗？

【参考方案3】：

# Mahalanobis Distance Matrix

import numpy as np
from scipy.spatial import distance
from scipy.spatial.distance import mahalanobis
from scipy.spatial.distance import pdist
from scipy.spatial.distance import squareform

# Example
Data = np.array([[1,2],[3,2],[4,3]])
Cov = np.cov(np.transpose([[1,2],[3,2],[4,3]]))
invCov = np.linalg.inv(Cov)

Y = pdist(Data, 'mahalanobis', invCov)
MD = squareform(Y) 

【讨论】：

您的答案可以通过额外的支持信息得到改进。请edit 添加更多详细信息，例如引用或文档，以便其他人可以确认您的答案是正确的。你可以找到更多关于如何写好答案的信息in the help center。