老饼讲解-机器学习 机器学习 神经网络 深度学习
相关基础数学

矩阵LU分解(三):Cholesky分解

作者 : 老饼 发表日期 : 2023-02-05 04:55:08 更新日期 : 2023-02-07 08:53:43
本站原创文章,转载请说明来自《老饼讲解-机器学习》www.bbbdata.com



本文先介绍什么是Cholesky矩阵分解及分解的方法

然后讲解相关公式的推导及及代码实现



    01. Cholesky矩阵分解的定义和分解流程    



本节讲解什么是Cholesky矩阵分解和如何对矩阵进行Cholesky分解



     Cholesky分解的定义    


Cholesky分解就是将正定对称矩阵分解成两个三角
A是正定对称矩阵,
求一个下三角矩阵 L,使得
 
即将A分解成

  , 




      Cholesky矩阵分解流程     


L的求解流程,是逐列求解,
即先求第1列,再求第2列,再求第3列....
 每列的求解分为两步:
  👉 先求对角元素    
  👉 再求非对角元素
 
 
假设前i-1列已求出,求第i列

求解 L第i列对角元素

 的求解公式如下

 
 即第i个对角元素为:根号(Aii-L第i行前(i-1)个元素的平方和)

 举例:
 

求解 L第i列所有元素

在求得后,第 i 列的求解公式如下
 
 


即L第i列为:(A的第i列- L前(i-1)列*各自第i行的元素)/
 


矩阵形式为:
  
 
✍️关于首列的计算公式
 特别地,第1列的计算公式为
   
 






      02. Cholesky矩阵分解流程原理推导    

 


本节讲解上节Cholesky分解流程中的公式是如何推导出来的



      Cholesky对角元素的公式推导     


对于L的对角元素,
公式的推导流程如下
             Aii等于L的 i 行乘的 i 列          
          的 i 列就是L的第i行的转置      
      i之后的元素是0,所以只要累加到i即可
 
 根据上式,即可得到

 






      Cholesky非角元素的公式推导     


假设前(i-1)列已求出,且也已求得,
现求L第i列的所有元素

           A的i列为 L乘以的 i 列                         
               L与互为转置                            
    A的i列为L每列乘以各自第i行元素后再相加    
 
   L的i行i列之后的元素为0,故不需i之后的列   
              独立抽出第i列  

 根据上式,即可得到
 
 


将上式写成矩阵形式,则为
  
 






      03. Cholesky矩阵分解的代码实现    


本节展示如何用用python代码实现Cholesky分解



    python代码实现Cholesky分解    


# -*- coding: utf-8 -*-
"""
 Cholesky分解的代码实现
"""
import numpy as np


def CholeskyDecompose(A):
    r,c=A.shape
    
    L=np.zeros((r,r))
    L[0,0] = np.sqrt(A[0,0])
    L[:,0] = A[:,0]/L[0,0]
    i=4
    for i in range(1,r):
        L[i,i] = np.sqrt(A[i,i]-sum(L[i,:i+1]*L[i,:i+1]))
        L[:,i] = (A[:,i]-L[:,:i]@L[i,:i].T)/L[i,i]
    return L

# 测试样例
if __name__ == "__main__":
    # 生成一个下三角矩阵L
    np.random.seed(0)
    rnd_mat = np.random.randint(0,10,(5,5))
    L_true = np.tril(rnd_mat,0)
    
    #将L*L^T组合成A
    A = L_true@L_true.T
    
    # 调用CholeskyDecompose对A进行Cholesky分解
    L = CholeskyDecompose(A)
    
    #打印结果
    print('----生成的L-----')
    print(L)
    print('\n----LL^T组合成的A-----')
    print(A)
    print('\n----从A分解得到的L-------')
    print(L)
代码运行结果如下

----生成的L-----
[[5. 0. 0. 0. 0.]
 [9. 3. 0. 0. 0.]
 [7. 6. 8. 0. 0.]
 [6. 7. 7. 8. 0.]
 [5. 9. 8. 9. 4.]]

----LL^T组合成的A-----
[[ 25  45  35  30  25]     
 [ 45  90  81  75  72]    
 [ 35  81 149 140 153] 
 [ 30  75 140 198 221] 
 [ 25  72 153 221 267]]

----从A分解得到的L-------
[[5. 0. 0. 0. 0.]
 [9. 3. 0. 0. 0.]
 [7. 6. 8. 0. 0.]
 [6. 7. 7. 8. 0.]
 [5. 9. 8. 9. 4.]]




    参考   


【1】Cholesky分解法 : https://blog.csdn.net/acdreamers/article/details/44656847










 End 






联系老饼