我目前正在研究内核方法，在某些时候我需要将一个非正半定矩阵(即相似性矩阵)制作成一个 PSD 矩阵。 我尝试了这种方法：
def makePSD(mat):
#make symmetric
k = (mat+mat.T)/2
#make PSD
min_eig = np.min(np.real(linalg.eigvals(mat)))
e = np.max([0, -min_eig + 1e-4])
mat = k + e*np.eye(mat.shape[0]);
return mat
但是如果我使用以下函数测试结果矩阵，它会失败：
def isPSD(A, tol=1e-8):
E,V = linalg.eigh(A)
return np.all(E >= -tol)
我还尝试了其他相关问题中建议的方法(如何计算最接近的正半定矩阵？)，但生成的矩阵也未能通过 isPSD 测试。
您对如何正确进行这种转换有什么建议吗？
我要说的第一件事是不要使用eigh来测试正确定性，因为eigh假设输入是埃尔米特的。这可能就是为什么你认为你引用的答案不起作用的原因。
我不喜欢这个答案，因为它有一个迭代(而且，我无法理解它的例子)，也不喜欢那里的另一个答案，它不承诺给你最好的正定矩阵，即最接近输入的弗罗贝尼乌斯范数(元素的平方和)。(我完全不知道你的问题中的代码应该做什么。
我确实喜欢Higham1988年论文的Matlab实现：https://www.mathworks.com/matlabcentral/fileexchange/42885-nearestspd 所以我把它移植到Python(针对Python 3进行了编辑更新)：
from numpy import linalg as la
import numpy as np
def nearestPD(A):
"""Find the nearest positive-definite matrix to input
A Python/Numpy port of John D'Errico's `nearestSPD` MATLAB code [1], which
credits [2].
[1] https://www.mathworks.com/matlabcentral/fileexchange/42885-nearestspd
[2] N.J. Higham, "Computing a nearest symmetric positive semidefinite
matrix" (1988): https://doi.org/10.1016/0024-3795(88)90223-6
"""
B = (A + A.T) / 2
_, s, V = la.svd(B)
H = np.dot(V.T, np.dot(np.diag(s), V))
A2 = (B + H) / 2
A3 = (A2 + A2.T) / 2
if isPD(A3):
return A3
spacing = np.spacing(la.norm(A))
# The above is different from [1]. It appears that MATLAB's `chol` Cholesky
# decomposition will accept matrixes with exactly 0-eigenvalue, whereas
# Numpy's will not. So where [1] uses `eps(mineig)` (where `eps` is Matlab
# for `np.spacing`), we use the above definition. CAVEAT: our `spacing`
# will be much larger than [1]'s `eps(mineig)`, since `mineig` is usually on
# the order of 1e-16, and `eps(1e-16)` is on the order of 1e-34, whereas
# `spacing` will, for Gaussian random matrixes of small dimension, be on
# othe order of 1e-16. In practice, both ways converge, as the unit test
# below suggests.
I = np.eye(A.shape[0])
k = 1
while not isPD(A3):
mineig = np.min(np.real(la.eigvals(A3)))
A3 += I * (-mineig * k**2 + spacing)
k += 1
return A3
def isPD(B):
"""Returns true when input is positive-definite, via Cholesky"""
try:
_ = la.cholesky(B)
return True
except la.LinAlgError:
return False
if __name__ == '__main__':
import numpy as np
for i in range(10):
for j in range(2, 100):
A = np.random.randn(j, j)
B = nearestPD(A)
assert (isPD(B))
print('unit test passed!')
除了找到最接近的正定矩阵外，上面的库还包括isPD，它使用 Cholesky 分解来确定矩阵是否为正定矩阵。这样，你不需要任何公差——任何想要正定的函数都会在其上运行 Cholesky，所以这是确定正定性的绝对最佳方法。
它最后还有一个基于蒙特卡洛的单元测试。如果您将其放入posdef.py并运行python posdef.py，它将在我的笔记本电脑上运行~一秒的单元测试。然后在代码中，您可以import posdef并调用posdef.nearestPD或posdef.isPD。
如果您这样做，代码也在 Gist 中。
我知道这个线程有点旧，但只是想说@user1231818链接的问题现在有一个令人满意的答案，至少在我测试过的情况下：https://stackoverflow.com/a/63131250/4733085
我将代码留在这里，但有关更多详细信息，请点击链接：
import numpy as np
def get_near_psd(A):
C = (A + A.T)/2
eigval, eigvec = np.linalg.eig(C)
eigval[eigval < 0] = 0
return eigvec.dot(np.diag(eigval)).dot(eigvec.T)
以防其他人面临同样的问题。 @tjiagoM和@Ahmed Fasih解释的上述两种方法是等价的，并且都应该给出相同的psd矩阵，以最小化Frobenius范数。 在论文的方程(2.1)和(2.2)中
'''
[2] N.J. Higham, "Computing a nearest symmetric positive semidefinite
matrix" (1988): https://doi.org/10.1016/0024-3795(88)90223-6
'''
在Ahmed的回答中，描述了tjiagoM使用的计算方法，实际上被称为计算最接近的psd矩阵的首选方法。