Lapack: Cholesky矩阵分解问题

Question

问题1

有人能推荐一种在python中进行Cholesky分解的不那么尴尬的方法吗？尤其是最后一行让我烦恼。

SigmaSqrt = matrix(Sigma)
cvxopt.lapack.potrf(SigmaSqrt)
SigmaSqrt = matrix(np.tril(SigmaSqrt))

问题2

我有一个问题，就是整行一列(例如，第一行中的所有元素和第一列中的所有元素)都是零，而失败的原因是矩阵不是正定的。处理这件事的最好方法是什么？

目前我正在做这件事：(这看起来很尴尬.)

try:
    SigmaSqrt = matrix(Sigma)
    cvxopt.lapack.potrf(SigmaSqrt)
    SigmaSqrt = matrix(np.tril(SigmaSqrt))
except ArithmeticError:
    SigmaSqrt = matrix(Sigma.ix[1:,1:])
    cvxopt.lapack.potrf(SigmaSqrt)
    SigmaSqrt = matrix(np.tril(SigmaSqrt))
    SigmaSqrt = sparse([[v0],[v0[1:].T, SigmaSqrt]])

Answer 1

你可以用numpy.linalg.cholesky。另外，如果一列或一行全部为零，则矩阵将是奇异的，至少在特征值上为零，因此不是正定的。由于Cholesky只定义为"Hermitian (对称的，如果实值的话)和正定“的矩阵，所以它不适用于它。

编辑：要“处理”您的问题取决于您想要什么。你做的任何事情都会使它发挥作用，这将是一个cholesky，它将不是原始矩阵的Cholesky。如果你正在进行一个迭代过程，并且它可以是模糊的，如果它已经是对称的，那么使用numpy.linalg.eigvalsh来寻找矩阵的特征值。设d是最负的特征值。然后设置A += (abs(d) + 1e-4) * np.indentity(len(A))。这将使它成为正定。

编辑：，这是Levenberg-Marquardt算法中使用的一个技巧。这链接到维基百科在牛顿法上的一篇文章，其中提到它是关于Levenberg-Marquard的文章没有进入this.Also，这里也是一篇关于它的文章。从根本上说，这将使(abs(d) + 1e-4)变换所有的特征值，从而使它们都为正，这是矩阵为正定的充分条件。

Answer 2

另一个选项是chompack模块：chompack主页 chompack.cholesky，我自己使用它与cvxopt模块结合使用。它与cvxopt的(稀疏)矩阵完美地工作。