求Python中半正定矩阵Cholesky分解的内置实现函数

NumPy和SciPy的原生cholesky函数只支持正定矩阵，处理半正定矩阵时会因奇异矩阵报错。以下是几种可行的解决思路：

• 矩阵平方根替代方案：scipy.linalg.sqrtm可直接计算半正定矩阵的平方根，结果L满足L @ L.T = A，效果等价于半正定矩阵的Cholesky分解。示例：

  import numpy as np
  from scipy.linalg import sqrtm
  
  # 构造半正定矩阵
  semi_pos_def = np.array([[2, 1], [1, 1]])
  L = sqrtm(semi_pos_def)
  # 验证分解正确性
  print(np.allclose(L @ L.T, semi_pos_def))
  

• 扰动近似分解：若允许微小误差，给半正定矩阵添加极小的对角扰动项，将其转为正定矩阵后使用标准Cholesky分解。示例：

  import numpy as np
  
  semi_pos_def = np.array([[2, 1], [1, 1]])
  n = semi_pos_def.shape[0]
  # 添加微小正则化项
  regularized = semi_pos_def + 1e-9 * np.eye(n)
  L = np.linalg.cholesky(regularized)
  print(np.allclose(L @ L.T, semi_pos_def, atol=1e-8))
  

• 注意事项：严格的Cholesky分解仅适用于正定矩阵，半正定矩阵可能不存在唯一的下三角分解。上述方案中，矩阵平方根是理论上更严谨的半正定矩阵分解方式，而扰动法是工程上常用的近似手段。

内容的提问来源于stack exchange，提问作者Ohad