Cholesky分解、乔列斯基分解

一、简介

1.1 定理

Cholesky分解法 又叫 平方根法，是一种分解 正定Hermite矩阵 (即 A = A H \boldsymbol A = \boldsymbol A^\mathrm H A=AH ) 的方法，以下我用实数域 (即 对称正定阵 A = A T \boldsymbol A = \boldsymbol A^\mathrm T A=AT) 来说明。

定理：若矩阵 A ∈ R n × n \boldsymbol A \in \mathbb R^{n\times n} A∈Rn×n 正定，且 A = A T \boldsymbol A=\boldsymbol A^{\mathrm T} A=AT，则存在唯一下三角矩阵 L ∈ R n × n \boldsymbol L \in \mathbb R^{n\times n} L∈Rn×n ，使得：
A = L L T \boldsymbol A = \boldsymbol L\boldsymbol L^{\mathrm T} A=LLT

证明：（暂略，有时间补）

1.2 分解方法

记:
A = ( a 11 a 12 … a 1 n a 21 a 22 … a 2 n ⋮ ⋮ ⋱ ⋮ a n 1 a n 2 … a n n )   L = ( l 11 l 21 l 22 ⋮ ⋮ ⋱ l n 1 l n 2 … l n n )    L T = ( l 11 l 21 … l n 1 l 22 … l n 2 ⋱ ⋮ l n n ) \boldsymbol A= \begin{pmatrix} a_{11} & a_{12} & \dots & a_{1n} \\ a_{21} & a_{22} & \dots & a_{2n} \\ \vdots & \vdots & \ddots & \vdots \\ a_{n1} & a_{n2} & \dots & a_{nn} \\ \end{pmatrix} \\ ~ \\ \boldsymbol L= \begin{pmatrix} l_{11} & \\ l_{21} & l_{22} \\ \vdots & \vdots & \ddots \\ l_{n1} & l_{n2} & \dots & l_{nn} \\ \end{pmatrix} ~~ \boldsymbol L^\mathrm T= \begin{pmatrix} l_{11} & l_{21} & \dots & l_{n1} \\ & l_{22} & \dots & l_{n2} \\ & & \ddots & \vdots \\ & & & l_{nn} \\ \end{pmatrix} A=⎝⎜⎜⎜⎛​a11​a21​⋮an1​​a12​a22​⋮an2​​……⋱…​a1n​a2n​⋮ann​​⎠⎟⎟⎟⎞​ L=⎝⎜⎜⎜⎛​l11​l21​⋮ln1​​l22​⋮ln2​​⋱…​lnn​​⎠⎟⎟⎟⎞​  LT=⎝⎜⎜⎜⎛​l11​​l21​l22​​……⋱​ln1​ln2​⋮lnn​​⎠⎟⎟⎟⎞​

那么：
( a 11 a 12 … a 1 n a 21 a 22 … a 2 n ⋮ ⋮ ⋱ ⋮ a n 1 a n 2 … a n n ) = ( l 11 l 21 l 22 ⋮ ⋮ ⋱ l n 1 l n 2 … l n n ) ( l 11 l 21 … l n 1 l 22 … l n 2 ⋱ ⋮ l n n )   = ( l 11 2 … … … … l 11 l 21 l 21 2 + l 22 2 … … … l 11 l 31 l 31 l 21 + l 32 l 22 l 31 2 + l 32 2 + l 33 2 … … ⋮ ⋮ ⋮ ⋱ ⋮ l 11 l n 1 l n 1 l 21 + l n 2 l 22 l n 1 l 31 + l n 2 l 32 + l n 3 l 33 … ∑ i = 1 n l n i 2 ) \begin{pmatrix} a_{11} & a_{12} & \dots & a_{1n} \\ a_{21} & a_{22} & \dots & a_{2n} \\ \vdots & \vdots & \ddots & \vdots \\ a_{n1} & a_{n2} & \dots & a_{nn} \\ \end{pmatrix}= \begin{pmatrix} l_{11} & \\ l_{21} & l_{22} \\ \vdots & \vdots & \ddots \\ l_{n1} & l_{n2} & \dots & l_{nn} \\ \end{pmatrix} \begin{pmatrix} l_{11} & l_{21} & \dots & l_{n1} \\ & l_{22} & \dots & l_{n2} \\ & & \ddots & \vdots \\ & & & l_{nn} \\ \end{pmatrix} \\ ~ \\ =\begin{pmatrix} l_{11}^2 & \dots & \dots & \dots & \dots \\ l_{11}l_{21} & l_{21}^2+l_{22}^2 & \dots & \dots & \dots \\ l_{11}l_{31} & l_{31}l_{21}+l_{32}l_{22} & l_{31}^2+l_{32}^2+l_{33}^2 & \dots & \dots \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ l_{11}l_{n1} & l_{n1}l_{21}+l_{n2}l_{22} & l_{n1}l_{31}+l_{n2}l_{32}+l_{n3}l_{33} & \dots & \displaystyle\sum_{i=1}^nl_{ni}^2 \\ \end{pmatrix} ⎝⎜⎜⎜⎛​a11​a21​⋮an1​​a12​a22​⋮an2​​……⋱…​a1n​a2n​⋮ann​​⎠⎟⎟⎟⎞​=⎝⎜⎜⎜⎛​l11​l21​⋮ln1​​l22​⋮ln2​​⋱…​lnn​​⎠⎟⎟⎟⎞​⎝⎜⎜⎜⎛​l11​​l21​l22​​……⋱​ln1​ln2​⋮lnn​​⎠⎟⎟⎟⎞​ =⎝⎜⎜⎜⎜⎜⎜⎜⎜⎛​l112​l11​l21​l11​l31​⋮l11​ln1​​…l212​+l222​l31​l21​+l32​l22​⋮ln1​l21​+ln2​l22​​……l312​+l322​+l332​⋮ln1​l31​+ln2​l32​+ln3​l33​​………⋱…​………⋮i=1∑n​lni2​​⎠⎟⎟⎟⎟⎟⎟⎟⎟⎞​

由于 A \boldsymbol A A 是对称矩阵，所以我们只看下三角：
( a 11 a 21 a 22 a 31 a 32 a 33 ⋮ ⋮ ⋮ ⋱ a n 1 a n 2 a n 3 … a n n ) = ( l 11 2 l 11 l 21 l 21 2 + l 22 2 l 11 l 31 l 31 l 21 + l 32 l 22 l 31 2 + l 32 2 + l 33 2 ⋮ ⋮ ⋮ ⋱ l 11 l n 1 l n 1 l 21 + l n 2 l 22 l n 1 l 31 + l n 2 l 32 + l n 3 l 33 … ∑ i = 1 n l n i 2 ) \begin{pmatrix} \textcolor{#FF0000}{a_{11}} \\ \textcolor{#FF0000}{a_{21}} & \textcolor{#DFCC00}{a_{22}} \\ \textcolor{#FF0000}{a_{31}} & \textcolor{#DFCC00}{a_{32}} & \textcolor{#00CC00}{a_{33}} \\ \textcolor{#FF0000}{\vdots} & \textcolor{#DFCC00}{\vdots} & \textcolor{#00CC00}{\vdots} & \textcolor{#0099FF}{\ddots} \\ \textcolor{#FF0000}{a_{n1}} & \textcolor{#DFCC00}{a_{n2}} & \textcolor{#00CC00}{a_{n3}} & \textcolor{#0099FF}{\dots} & \textcolor{#BB00FF}{a_{nn}} \\ \end{pmatrix}= \begin{pmatrix} \textcolor{#FF0000}{l_{11}^2} \\ \textcolor{#FF0000}{l_{11}l_{21}} & \textcolor{#DFCC00}{l_{21}^2+l_{22}^2} \\ \textcolor{#FF0000}{l_{11}l_{31}} & \textcolor{#DFCC00}{l_{31}l_{21}+l_{32}l_{22}} & \textcolor{#00CC00}{l_{31}^2+l_{32}^2+l_{33}^2} \\ \textcolor{#FF0000}{\vdots} & \textcolor{#DFCC00}{\vdots} & \textcolor{#00CC00}{\vdots} & \textcolor{#0099FF}{\ddots} \\ \textcolor{#FF0000}{l_{11}l_{n1}} & \textcolor{#DFCC00}{l_{n1}l_{21}+l_{n2}l_{22}} & \textcolor{#00CC00}{l_{n1}l_{31}+l_{n2}l_{32}+l_{n3}l_{33}} & \textcolor{#0099FF}{\dots} & \textcolor{#BB00FF}{\displaystyle\sum_{i=1}^nl_{ni}^2} \\ \end{pmatrix} ⎝⎜⎜⎜⎜⎜⎛​a11​a21​a31​⋮an1​​a22​a32​⋮an2​​a33​⋮an3​​⋱…​ann​​⎠⎟⎟⎟⎟⎟⎞​=⎝⎜⎜⎜⎜⎜⎜⎜⎜⎛​l112​l11​l21​l11​l31​⋮l11​ln1​​l212​+l222​l31​l21​+l32​l22​⋮ln1​l21​+ln2​l22​​l312​+l322​+l332​⋮ln1​l31​+ln2​l32​+ln3​l33​​⋱…​i=1∑n​lni2​​⎠⎟⎟⎟⎟⎟⎟⎟⎟⎞​

按照从左到右 (红黄绿蓝紫) 的顺序，逐列对应元素计算，便可将所有 L \boldsymbol L L 的元素 l i j l_{ij} lij​ 求出来

二、应用

2.1 线性方程求解

对于线性方程组 A X = B \boldsymbol A\boldsymbol X=\boldsymbol B AX=B，其中 A ∈ R n × n \boldsymbol A \in \mathbb R^{n\times n} A∈Rn×n 正定，且 A = A T \boldsymbol A=\boldsymbol A^{\mathrm T} A=AT，那么求解方程可以使用 Cholesky分解：

1. 对矩阵 A \boldsymbol A A 进行 Cholesky分解 得， A = L L T \boldsymbol A = \boldsymbol L\boldsymbol L^{\mathrm T} A=LLT，则原方程化为 L L T X = B \boldsymbol L\boldsymbol L^{\mathrm T}\boldsymbol X=\boldsymbol B LLTX=B
2. 令 L T X = Y \boldsymbol L^{\mathrm T}\boldsymbol X=\boldsymbol Y LTX=Y，此时，解下三角方程 L Y = B \boldsymbol L\boldsymbol Y=\boldsymbol B LY=B，求得 Y \boldsymbol Y Y
3. 解上三角方程 L T X = Y \boldsymbol L^{\mathrm T}\boldsymbol X=\boldsymbol Y LTX=Y，求得 X \boldsymbol X X

通过 Cholesky分解 将 普通线性方程求解 改为两次简单的 三角阵方程组求解 ，降低计算复杂度。

2.2 求逆矩阵

三角矩阵的逆比较好求，从而可以很快求出原矩阵的逆。

三、扩展Cholesky分解

3.1 简介

从 1.2 节可以知道，矩阵 L \boldsymbol L L 对角线上的元素在计算时，都需要开方操作，增加计算量的同时还可能损失精度。引进一种 Cholesky的扩展分解方法：
A = L D L T \boldsymbol A = \boldsymbol L\boldsymbol D\boldsymbol L^{\mathrm T} A=LDLT

证明：
我们已知，正定对称阵 A \boldsymbol A A 可以分解为： A = L L T \boldsymbol A = \boldsymbol L\boldsymbol L^{\mathrm T} A=LLT

而 L \boldsymbol L L 可以写成：

L = ( l 11 l 21 l 22 ⋮ ⋮ ⋱ l n 1 l n 2 … l n n )   = ( 1 l 21 / l 11 1 ⋮ ⋮ ⋱ l n 1 / l 11 l n 2 / l 22 … 1 ) ( l 11 l 22 ⋱ l n n )   = d e f L 1 Λ \begin{aligned} \boldsymbol L &=\begin{pmatrix} l_{11} & \\ l_{21} & l_{22} \\ \vdots & \vdots & \ddots \\ l_{n1} & l_{n2} & \dots & l_{nn} \\ \end{pmatrix} \\ ~ \\ &=\begin{pmatrix} 1 & \\ l_{21}/l_{11} & 1 \\ \vdots & \vdots & \ddots \\ l_{n1}/l_{11} & l_{n2}/l_{22} & \dots & 1 \\ \end{pmatrix} \begin{pmatrix} l_{11} & \\ & l_{22} \\ & & \ddots \\ & & & l_{nn} \\ \end{pmatrix} \\ ~ \\ &\overset{\mathrm{def}}{=}\boldsymbol L_1 \boldsymbol \Lambda \end{aligned} L  ​=⎝⎜⎜⎜⎛​l11​l21​⋮ln1​​l22​⋮ln2​​⋱…​lnn​​⎠⎟⎟⎟⎞​=⎝⎜⎜⎜⎛​1l21​/l11​⋮ln1​/l11​​1⋮ln2​/l22​​⋱…​1​⎠⎟⎟⎟⎞​⎝⎜⎜⎛​l11​​l22​​⋱​lnn​​⎠⎟⎟⎞​=defL1​Λ​

同理 L T = Λ L 1 T \boldsymbol L^\mathrm T = \boldsymbol \Lambda\boldsymbol L_1^\mathrm T LT=ΛL1T​

那么 A = L L T = L 1 Λ 2 L 1 T = L 1 D L 1 T \boldsymbol A = \boldsymbol L\boldsymbol L^{\mathrm T}=\boldsymbol L_1\boldsymbol \Lambda^2\boldsymbol L_1^{\mathrm T} = \boldsymbol L_1\boldsymbol D\boldsymbol L_1^{\mathrm T} A=LLT=L1​Λ2L1T​=L1​DL1T​

3.2 分解方法

记:
A = ( a 11 a 12 … a 1 n a 21 a 22 … a 2 n ⋮ ⋮ ⋱ ⋮ a n 1 a n 2 … a n n )   L = ( 1 l 21 1 ⋮ ⋮ ⋱ l n 1 l n 2 … 1 )    L T = ( 1 l 21 … l n 1 1 … l n 2 ⋱ ⋮ 1 )    D = ( d 1 d 2 ⋱ d n ) \boldsymbol A= \begin{pmatrix} a_{11} & a_{12} & \dots & a_{1n} \\ a_{21} & a_{22} & \dots & a_{2n} \\ \vdots & \vdots & \ddots & \vdots \\ a_{n1} & a_{n2} & \dots & a_{nn} \\ \end{pmatrix} \\ ~ \\ \boldsymbol L= \begin{pmatrix} 1 & \\ l_{21} & 1 \\ \vdots & \vdots & \ddots \\ l_{n1} & l_{n2} & \dots & 1 \\ \end{pmatrix} ~~ \boldsymbol L^\mathrm T= \begin{pmatrix} 1 & l_{21} & \dots & l_{n1} \\ & 1 & \dots & l_{n2} \\ & & \ddots & \vdots \\ & & & 1 \\ \end{pmatrix} ~~ \boldsymbol D= \begin{pmatrix} d_1 & & & \\ & d_2 & & \\ & & \ddots & \\ & & & d_n \\ \end{pmatrix} A=⎝⎜⎜⎜⎛​a11​a21​⋮an1​​a12​a22​⋮an2​​……⋱…​a1n​a2n​⋮ann​​⎠⎟⎟⎟⎞​ L=⎝⎜⎜⎜⎛​1l21​⋮ln1​​1⋮ln2​​⋱…​1​⎠⎟⎟⎟⎞​  LT=⎝⎜⎜⎜⎛​1​l21​1​……⋱​ln1​ln2​⋮1​⎠⎟⎟⎟⎞​  D=⎝⎜⎜⎛​d1​​d2​​⋱​dn​​⎠⎟⎟⎞​

那么：
( a 11 a 12 … a 1 n a 21 a 22 … a 2 n ⋮ ⋮ ⋱ ⋮ a n 1 a n 2 … a n n ) = ( 1 l 21 1 ⋮ ⋮ ⋱ l n 1 l n 2 … 1 ) ( d 1 d 2 ⋱ d n ) ( 1 l 21 … l n 1 1 … l n 2 ⋱ ⋮ 1 )   = ( d 1 … … … … d 1 l 21 d 1 l 21 2 + d 2 … … … d 1 l 31 d 1 l 31 l 21 + l 32 l 22 d 1 l 31 2 + d 2 l 32 2 + d 3 … … ⋮ ⋮ ⋮ ⋱ ⋮ d 1 l n 1 d 1 l n 1 l 21 + d 2 l n 2 d 1 l n 1 l 31 + d 2 l n 2 l 32 + d 3 l n 3 … d n + ∑ i = 1 n − 1 d i l n i 2 ) \begin{pmatrix} a_{11} & a_{12} & \dots & a_{1n} \\ a_{21} & a_{22} & \dots & a_{2n} \\ \vdots & \vdots & \ddots & \vdots \\ a_{n1} & a_{n2} & \dots & a_{nn} \\ \end{pmatrix}= \begin{pmatrix} 1 & \\ l_{21} & 1 \\ \vdots & \vdots & \ddots \\ l_{n1} & l_{n2} & \dots & 1 \\ \end{pmatrix} \begin{pmatrix} d_1 & & & \\ & d_2 & & \\ & & \ddots & \\ & & & d_n \\ \end{pmatrix} \begin{pmatrix} 1 & l_{21} & \dots & l_{n1} \\ & 1 & \dots & l_{n2} \\ & & \ddots & \vdots \\ & & & 1 \\ \end{pmatrix} \\ ~ \\ =\begin{pmatrix} d_1 & \dots & \dots & \dots & \dots \\ d_1l_{21} & d_1l_{21}^2+d_2 & \dots & \dots & \dots \\ d_1l_{31} & d_1l_{31}l_{21}+l_{32}l_{22} & d_1l_{31}^2+d_2l_{32}^2+d_3 & \dots & \dots \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ d_1l_{n1} & d_1l_{n1}l_{21}+d_2l_{n2} & d_1l_{n1}l_{31}+d_2l_{n2}l_{32}+d_3l_{n3} & \dots & \displaystyle d_n+\sum_{i=1}^{n-1}d_il_{ni}^2 \\ \end{pmatrix} ⎝⎜⎜⎜⎛​a11​a21​⋮an1​​a12​a22​⋮an2​​……⋱…​a1n​a2n​⋮ann​​⎠⎟⎟⎟⎞​=⎝⎜⎜⎜⎛​1l21​⋮ln1​​1⋮ln2​​⋱…​1​⎠⎟⎟⎟⎞​⎝⎜⎜⎛​d1​​d2​​⋱​dn​​⎠⎟⎟⎞​⎝⎜⎜⎜⎛​1​l21​1​……⋱​ln1​ln2​⋮1​⎠⎟⎟⎟⎞​ =⎝⎜⎜⎜⎜⎜⎜⎜⎜⎛​d1​d1​l21​d1​l31​⋮d1​ln1​​…d1​l212​+d2​d1​l31​l21​+l32​l22​⋮d1​ln1​l21​+d2​ln2​​……d1​l312​+d2​l322​+d3​⋮d1​ln1​l31​+d2​ln2​l32​+d3​ln3​​………⋱…​………⋮dn​+i=1∑n−1​di​lni2​​⎠⎟⎟⎟⎟⎟⎟⎟⎟⎞​

计算的方式与之前类似，从左至右直接逐列对应计算即可。

可以看出来这种方式不需要开方。