在尝试对周期边界条件的2D数组的方差-协方差矩阵执行cholesky分解时,在某些参数组合下,我总是得到LinAlgError: Matrix is not positive definite - Cholesky decomposition cannot be computed
。由于脚本很简单,因此不确定是numpy.linalg
还是实现问题:
sigma = 3.
U = 4
def FromListToGrid(l_):
i = np.floor(l_/U)
j = l_ - i*U
return np.array((i,j))
Ulist = range(U**2)
Cov = []
for l in Ulist:
di = np.array([np.abs(FromListToGrid(l)[0]-FromListToGrid(i)[0]) for i, x in enumerate(Ulist)])
di = np.minimum(di, U-di)
dj = np.array([np.abs(FromListToGrid(l)[1]-FromListToGrid(i)[1]) for i, x in enumerate(Ulist)])
dj = np.minimum(dj, U-dj)
d = np.sqrt(di**2+dj**2)
Cov.append(np.exp(-d/sigma))
Cov = np.vstack(Cov)
W = np.linalg.cholesky(Cov)
尝试删除潜在的奇点也未能解决问题。任何帮助深表感谢。
最佳答案
深入研究问题,我尝试打印Cov矩阵的特征值。
print np.linalg.eigvalsh(Cov)
答案是这样的
[-0.0801339 -0.0801339 0.12653595 0.12653595 0.12653595 0.12653595 0.14847999 0.36269785 0.36269785 0.36269785 0.36269785 1.09439988 1.09439988 1.09439988 1.09439988 9.6772531 ]
啊哈!注意前两个负特征值吗?现在,当且仅当其所有特征值均为正时,矩阵才为正定。因此,矩阵的问题不在于它接近“零”,而在于它是“负”。为了扩展@duffymo类比,这是线性代数,等效于尝试取负数的平方根。
现在,让我们尝试执行相同的操作,但是这次是scipy。
scipy.linalg.cholesky(Cov, lower=True)
而且这还不能说更多
numpy.linalg.linalg.LinAlgError: 12-th leading minor not positive definite
这说明了更多信息(尽管我无法真正理解为什么它提示第12个未成年人)。
底线,矩阵不是很接近“零”,而是更像“负”