Choleského rozklad

Choleského rozklad matice $𝑨$ na součin dolní trojúhelníkové matice $𝑳$ a její transpozice.

Choleského rozklad (také Choleského dekompozice nebo Choleského faktorizace) je metoda lineární algebry, kterou lze každou reálnou pozitivně definitní matici rozložit na součin dolní trojúhelníkové matice a její transpozice. Obecněji lze pojem zavést i pro komplexní matice.

Rozklad je pojmenován po francouzském matematikovi André-Louisovi Choleském (1875–1918, výslovnost [šoleski], Šablona:IPA2), který jej vyvinul před rokem 1914 při triangulaci Kréty francouzskou Service géographique de l’armée.

Pro řešení soustav lineárních rovnic s pozitivně definitní maticí je Choleského rozklad zhruba dvakrát efektivnější než LU rozklad.

Definice

Pro každou pozitivně semidefinitní komplexní matici $𝑳$ existuje dolní trojúhelníková matice $𝑳$ taková, že platí:

$𝑨={𝑳𝑳}^{\mathrm{H}}$

Uvedený zápis matice $𝑨$ jako součin ${𝑳𝑳}^{\mathrm{H}}$ se nazývá Choleského rozklad matice $𝑨$. Dolní trojúhelníková matice $𝑳$ se nazývá Choleského faktor matice $𝑨$. Symbol ${𝑳}^{\mathrm{H}}$ značí matici hermitovsky sdruženou k matici $𝑳$ (též nazývanou hermitovská transpozice).

Ten samý rozklad lze zapsat ve tvaru $𝑨={𝑼}^{\mathrm{H}}𝑼$, kde $𝑼$ je horní trojúhelníková, neboť $𝑼={𝑳}^{\mathrm{H}}$.

Reálné pozitivně definitní matice mají Choleského faktory reálné. Pro ně platí: ${𝑳}^{\mathrm{H}}={𝑳}^{\mathrm{T}}$, a proto lze rozklad zapsat ve tvaru:

$𝑨={𝑳𝑳}^{\mathrm{T}}$

Ukázka

Symetrická reálná matice

$𝑨=\left(\begin{array}{ccc}4& 12& -16\\ 12& 37& -43\\ -16& -43& 98\end{array}\right)$

má Choleského rozklad ${𝑳𝑳}^{\mathrm{T}}$:

$\left(\begin{array}{ccc}4& 12& -16\\ 12& 37& -43\\ -16& -43& 98\end{array}\right)=\left(\begin{array}{ccc}2& 0& 0\\ 6& 1& 0\\ -8& 5& 3\end{array}\right)\left(\begin{array}{ccc}2& 6& -8\\ 0& 1& 5\\ 0& 0& 3\end{array}\right)$

Vlastnosti

• Choleského faktor $𝑳$ je regulární právě když je daná matice $𝑨$ regulární.

• Má-li matice $𝑨$ Choleského rozklad ${𝑳𝑳}^{\mathrm{H}}$, je hermitovská, resp. u reálných symetrická, protože ${𝑨}^{\mathrm{H}}=({𝑳𝑳}^{\mathrm{H}}{)}^{\mathrm{H}}=({𝑳}^{\mathrm{H}}{)}^{\mathrm{H}}{𝑳}^{\mathrm{H}}={𝑳𝑳}^{\mathrm{H}}=𝑨$.

• Má-li matice $𝑨$ Choleského rozklad s regulárním Choleského faktorem $𝑳$ je pozitivně definitní. Pro libovolné $𝒙\ne 0$ vyplývá z regularity matice $𝑳$, že také ${𝑳}^{\mathrm{H}}𝒙\ne 0$, a potom

${𝒙}^{\mathrm{H}}𝑨𝒙={𝒙}^{\mathrm{H}}{𝑳𝑳}^{\mathrm{H}}𝒙=({𝑳}^{\mathrm{H}}𝒙{)}^{\mathrm{H}}{𝑳}^{\mathrm{H}}𝒙=⟨{𝑳}^{\mathrm{H}}𝒙|{𝑳}^{\mathrm{H}}𝒙⟩>0$, přičemž v předposlední výraz $⟨\cdot |\cdot ⟩$ značí standardní skalární součin na ${ℂ}^n$.

• Choleského rozklad není jednoznačný, např. matici $\left(\begin{array}{cc}1& 1\\ 1& 2\end{array}\right)$ lze rozložit čtyřmi způsoby s Choleského faktory: $\left(\begin{array}{cc}1& 0\\ 1& 1\end{array}\right),\left(\begin{array}{cc}1& 0\\ 1& -1\end{array}\right),\left(\begin{array}{cc}-1& 0\\ -1& 1\end{array}\right)$ a $\left(\begin{array}{cc}-1& 0\\ -1& -1\end{array}\right)$.

• Choleského faktory pozitivně semidefinitních (i komplexních) matic mají na diagonále vždy reálná čísla.

• Pouze jeden z Choleského faktorů pozitivně definitních matic má všechny prvky na diagonále kladné.

• Pokud je hermitovská matice $𝑨$ pouze pozitivně semidefinitní, a nikoli pozitivně definitní, pak má stále Choleského rozklad, kde alespoň jeden prvek na diagonále $𝑳$ je nulový. Choleského faktorů může být i nekonečně mnoho, například rozkladem matice $\left(\begin{array}{cc}0& 0\\ 0& 1\end{array}\right)$ je každá matice $\left(\begin{array}{cc}0& 0\\ \cos \phi & \sin \phi \end{array}\right)$, kde $\phi \in (0,\pi )$.

• Mezi Choleského faktory pozitivně semidefinitních matic hodnosti $r$ lze nalézt právě jeden takový, že má $r$ kladných prvků na diagonále a $n-r$ sloupců se samými nulami. Jinak řečeno, v tomto případě existuje alespoň jedna permutační matice $𝑷$ taková, že matice ${𝑷𝑨𝑷}^{\mathrm{T}}$ má jednoznačný Choleského rozklad ve tvaru $𝑳=\left(\begin{array}{cc}{𝑳}_1& 0\\ {𝑳}_2& 0\end{array}\right)$, kde ${𝑳}_1$ je dolní trojúhelníková matice hodnosti $r$ s kladnou diagonálou.

LDL rozklad

S Choleského rozkladem úzce souvisí rozklad dané matice na součin:

$𝑨={𝑳𝑫𝑳}^{\mathrm{H}}$,

kde $𝑳$ je dolní trojúhelníková s 1 na diagonále a $𝑫$ je diagonální.

LDL rozklad lze vypočítat a použít v podstatě stejnými algoritmy jako klasický Choleského rozklad, ovšem bez použití odmocnin.

Ukázka

Matice $𝑨$ z předchozí ukázky má LDL rozklad:

$\left(\begin{array}{ccc}4& 12& -16\\ 12& 37& -43\\ -16& -43& 98\end{array}\right)=\left(\begin{array}{ccc}1& 0& 0\\ 3& 1& 0\\ -4& 5& 1\end{array}\right)\left(\begin{array}{ccc}4& 0& 0\\ 0& 1& 0\\ 0& 0& 9\end{array}\right)\left(\begin{array}{ccc}1& 3& -4\\ 0& 1& 5\\ 0& 0& 1\end{array}\right)$

Choleského faktor z předchozí ukázky lze spočítat pomocí součinu s odmocninou z diagonální matice:

$\left(\begin{array}{ccc}2& 0& 0\\ 6& 1& 0\\ -8& 5& 3\end{array}\right)=\left(\begin{array}{ccc}1& 0& 0\\ 3& 1& 0\\ -4& 5& 1\end{array}\right)\left(\begin{array}{ccc}2& 0& 0\\ 0& 1& 0\\ 0& 0& 3\end{array}\right)$

LDL rozklad může mít například i matice, která je negativně semidefinitní.

$\left(\begin{array}{ccc}-2& -4& 2\\ -4& -9& 4\\ 2& 4& -2\end{array}\right)=\left(\begin{array}{ccc}1& 0& 0\\ 2& 1& 0\\ -1& 0& 1\end{array}\right)\left(\begin{array}{ccc}-2& 0& 0\\ 0& -1& 0\\ 0& 0& 0\end{array}\right)\left(\begin{array}{ccc}1& 2& -1\\ 0& 1& 0\\ 0& 0& 1\end{array}\right)$

Vlastnosti

• Je-li $𝑨$ pozitivně definitní, pak jsou všechny prvky na diagonále $𝑫$ kladné. Z LDL rozkladu lze pak odvodit klasický Choleského rozklad s faktorem $𝑪=𝑳{𝑫}^{1/2}$ pomocí vztahu:

$𝑨={𝑳𝑫𝑳}^{\mathrm{H}}=𝑳{𝑫}^{1/2}{\left({𝑫}^{1/2}\right)}^{\mathrm{H}}{𝑳}^{\mathrm{H}}=𝑳{𝑫}^{1/2}{\left(𝑳{𝑫}^{1/2}\right)}^{\mathrm{H}}=𝑪{𝑪}^{\mathrm{H}}$

• Naopak, má-li pozitivně definitní matice klasický Choleského rozklad $𝑨=𝑪{𝑪}^{\mathrm{H}}$, a matice $𝑺$ je diagonální matice, která obsahuje hlavní diagonálu $𝑪$, pak $𝑨$ lze rozložit jako $𝑳𝑫{𝑳}^{\mathrm{H}}$, kde:

$𝑳=𝑪{𝑺}^{-1}$, tím se sloupce naškálují tak, aby prvky na diagonále byly rovny 1,

$𝑫={𝑺}^2$.

• Pozitivně semidefinitní matice mají LDL rozklad právě když se hodnosti matic $𝑫$ a $𝑨$ shodují.

• Pro existenci LDL rozkladu hermitovské (ne nutně pozitivně definitní) matice například stačí, aby prvních $n-1$ hlavních vedoucích minorů matice $𝑨$ bylo nenulových.

• Není-li matice pozitivně semidefinitní, čili je negativně (semi)definitní nebo indefinitní, a přitom má LDL rozklad, potom se na diagonále $𝑫$ vyskytne alespoň jedno záporné číslo.

• Matice $𝑨$ a $𝑫$ mají shodný determinant a ten je roven součinu prvků na diagonále matice $𝑫$.

Výpočet

Z rozepsání součinu pro matice řádu 3

$\begin{array}{cc}𝑨={𝑳𝑳}^T& =\left(\begin{array}{ccc}{l}_{11}& 0& 0\\ l_{21}& l_{22}& 0\\ l_{31}& l_{32}& l_{33}\end{array}\right)\left(\begin{array}{ccc}{l}_{11}& l_{21}& l_{31}\\ 0& l_{22}& l_{32}\\ 0& 0& l_{33}\end{array}\right)\\ & =\left(\begin{array}{ccc}{l}_{11}^2& & (\text{symetricky})\\ l_{21}{l}_{11}& l_{21}^2+l_{22}^2& \\ l_{31}{l}_{11}& l_{31}{l}_{21}+l_{32}{L}_{22}& l_{31}^2+l_{32}^2+l_{33}^2\end{array}\right),\end{array}$

vyplývá, že Choleského faktor s kladnou diagonálou je dán výrazem:

$\begin{array}{c}𝑳=\left(\begin{array}{ccc}\sqrt{a_{11}}& 0& 0\\ a_{21}/l_{11}& \sqrt{a_{22}-l_{21}^2}& 0\\ a_{31}/l_{11}& (a_{32}-l_{31}{l}_{21})/l_{22}& \sqrt{a_{33}-l_{31}^2-l_{32}^2}\end{array}\right)\end{array}$

Obecně je možné prvky matice $𝑳$ počítat po sloupcích zleva doprava a v každém sloupci odshora dolů.

Pro první sloupec platí následující.

$a_{11}=l_{11}{l}_{11}⟶l_{11}=\sqrt{a_{11}}$

$a_{21}=l_{21}{l}_{11}⟶l_{21}=a_{21}/l_{11}$

$⋮$

$a_{n1}=l_{n1}{l}_{11}⟶l_{n1}=a_{n1}/l_{11}$

Pro druhý sloupec platí:

$a_{22}=l_{21}{l}_{21}+l_{22}{l}_{22}⟶l_{22}=\sqrt{a_{22}-l_{21}^2}$

$a_{32}=l_{31}{l}_{21}+l_{32}{l}_{22}⟶l_{32}=(a_{32}-l_{31}{l}_{21})/l_{22}$

$⋮$

$a_{n2}=l_{n1}{l}_{21}+l_{n2}{l}_{22}⟶l_{n2}=(a_{n2}-l_{n1}{l}_{21})/l_{22}$

Pro prvky na diagonále lze, vzhledem ke znalosti celého řádku vlevo od prvku, odvodit následující vzorec:

$a_{ii}=\sum _{k=1}^i{l}_{ik}^2⟶l_{ii}=(\pm )\sqrt{a_{ii}-\sum _{k=1}^{i-1}{l}_{ik}^2}$

Pro prvky pod diagonálou vyplývá podobně následující vztah:

$a_{ij}=\sum _{k=1}^j{l}_{ik}{l}_{jk}⟶l_{ij}=\frac{1}{l_{jj}}(a_{ij}-\sum _{k=1}^{j-1}{l}_{ik}{l}_{jk})$ pro $i>j$

U prvků na diagonále je možné vzít hodnotu odmocniny se záporným znaménkem, což vyvolá změnu na nich závisejících prvků mimo diagonálu.

Pro komplexní pozitivně definitní matice platí analogické vztahy (pruhem je značeno komplexně sdružené číslo):

$a_{ii}=\sum _{k=1}^i{l}_{ik}\overline{l_{ik}}⟶l_{ii}=(\pm )\sqrt{a_{ii}-\sum _{k=1}^{i-1}{l}_{ik}\overline{l_{ik}}}$

$a_{ij}=\sum _{k=1}^j{l}_{ik}\overline{l_{jk}}⟶l_{ij}=\frac{1}{l_{jj}}(a_{ij}-\sum _{k=1}^{j-1}{l}_{ik}\overline{l_{jk}})$ pro $i>j$

Výraz pod odmocninou je pro pozitivně definitní matice vždy kladný.

Vzor čtení (bíle) a zápisu (žlutě) pro výpočet Choleského rozkladu na místě podle Choleského–Banachiewiczova algoritmu pro matici řádu 5

Pseudokód

Výpočty ve výše uvedených vzorcích lze provádět různými způsoby. Varianta pojmenovaná po Tadeuszi Banachiewiczovi vypočte dolní trojúhelníkovou matici řádek po řádku a přitom na místě. V pseudokódu je uveden postup rozkladu matice $𝑨$ do tvaru ${𝑳𝑳}^{\mathrm{H}}$:

  Input: hermitovská matice A řádu n reprezentovaná svou dolní trojúhelníkovou polovinou
  Output: dolní trojúhelníková část Choleského faktoru L 
  For i = 1 To n
    For j = 1 To i
      Suma = a(i, j)
      For k = 1 To j-1
        Suma = Suma - a(i, k) * conj(a(j, k))
      If i > j Then
        a(i, j) = Suma / a(j, j)  // Prvek je pod diagonálou.
      Else If Suma > 0 Then       // Prvek na diagonále
        a(i, i) = Sqrt(Suma)      // … musí být vždy nezáporný.
      Else
        ERROR("Matice není pozitivně definitní.")
  Return: L=A

Algoritmus pracuje na místě: postupně mění matici $𝑨$ na $𝑳$, aniž by bylo třeba alokovat další paměť pro zápis výsledné matice. Navíc využívá pouze dolní trojúhelníkovou matici, protože hodnoty prvků nad diagonálou lze dopočítat s využitím vlastnosti, že daná matice $𝑨$ je hermitovská. Výsledný Choleského faktor $𝑳$ je třeba vzít tak, že má prvky nad diagonálou nulové.

Výpočetní složitost

Časová složitost běžně používaných algoritmů pro výpočet Choleského rozkladu je obecně $O(n^3)$.  Přesněji, na reálných maticích jde o $\frac{1}{3}{n}^3+O(n^2)$ aritmetických operací s prvky dané matice, konkrétně $\frac{n^3}{6}$ součinů i součtů, $\frac{n^2}{2}$ dělení a $n$ odmocnin. Komplexní matice oproti tomu vyžadují $\frac{4}{6}{n}^3$ součinů i součtů.

Pro srovnání, LU rozklad coby implementace Gaussovy eliminace, vyžaduje přibližně dvakrát více aritmetických operací.

Numerické záležitosti

Choleského rozklad je bezpodmínečně zpětně stabilní.

Je-li daná matice pozitivně definitní, jsou čísla pod odmocninami vždy kladná v přesné aritmetice. Zaokrouhlovací chyby mohou tuto vlastnost porušit a v takovém případě algoritmus nemůže pokračovat. Tento případ však může nastat, jen je-li matice velmi špatně podmíněna.

LDL rozklad

Výpočtu odmocnin se lze vyhnout ve výpočtu LDL rozkladu. Ten lze spočítat i v přesné zlomkové aritmetice, jak lze odvodit následovně. Pro rozklad reálné matice řádu 3 platí:

$\begin{array}{cc}𝑨={𝑳𝑫𝑳}^{\mathrm{T}}& =\left(\begin{array}{ccc}1& 0& 0\\ l_{21}& 1& 0\\ l_{31}& l_{32}& 1\end{array}\right)\left(\begin{array}{ccc}{d}_{11}& 0& 0\\ 0& d_{22}& 0\\ 0& 0& d_{33}\end{array}\right)\left(\begin{array}{ccc}1& l_{21}& l_{31}\\ 0& 1& l_{32}\\ 0& 0& 1\end{array}\right)\\ & =\left(\begin{array}{ccc}{d}_{11}& & (\mathrm{s}\mathrm{y}\mathrm{m}\mathrm{e}\mathrm{t}\mathrm{r}\mathrm{i}\mathrm{c}\mathrm{k}\mathrm{y})\\ l_{21}{d}_{11}& l_{21}^2{d}_{11}+d_{22}& \\ l_{31}{d}_{11}& l_{31}{l}_{21}{d}_{11}+l_{32}{d}_{22}& l_{31}^2{d}_{11}+l_{32}^2{d}_{22}+d_{33}.\end{array}\right).\end{array}$

Obecně jsou prvky matic $𝑫$ a $𝑳$ i vyšších řádů dány následujícími rekurentními vzorci:

$d_{jj}=a_{jj}-\sum _{k=1}^{j-1}{l}_{jk}^2{d}_{kk},$

$l_{ij}=\frac{1}{d_{jj}}(a_{ij}-\sum _{k=1}^{j-1}{l}_{ik}{l}_{jk}{d}_{kk})$ pro $i>j$.

Pro komplexní matice je třeba výrazy na pravé straně upravit následovně:

$d_{jj}=A_{jj}-\sum _{k=1}^{j-1}{l}_{jk}\overline{l_{jk}}{d}_{kk},$

$l_{ij}=\frac{1}{d_{jj}}(a_{ij}-\sum _{k=1}^{j-1}{l}_{ik}\overline{l_{jk}}{d}_{kk})$ pro $i>j$.

Vzorec přístupu k prvkům matice opět umožňuje, aby byl v případě potřeby celý výpočet proveden na místě.

Aplikace

Numerické řešení soustavy lineárních rovnic

Choleského rozklad se používá především pro numerické řešení lineárních rovnic $𝑨𝒙=𝒃$ s pozitivně definitní maticí soustavy a to tak, že se nejprve provede Choleského rozklad $𝑨={𝑳𝑳}^{\mathrm{H}}$, potom se dopřednou substitucí určí řešení $𝒚$ soustavy $𝑳𝒚=𝒃$ a nakonec se zpětnou substitucí vyřeší soustava ${𝑳}^{\mathrm{H}}𝒙=𝒚$.

Vzhledem k tomu, že matice obou soustav jsou trojúhelníkové, je řešení uvedených soustav snadné. Choleského rozklad (nebo jeho LDL varianta, kde ani není třeba odmocňovat) je u těchto soustav oblíbenou pro svou účinnost a numerickou stabilitu. Ve srovnání s LU rozkladem je zhruba dvakrát efektivnější.

Inverzní matice

Matici inverzní k pozitivně definitní matici lze spočítat pomocí Choleského rozkladu podobným způsobem jako při řešení soustav lineárních rovnic v čase $O(n^3)$. Postup lze provést i na místě.

Libovolná komplexní regulární matice $𝑩$ může být invertována pomocí následující identity, protože ${𝑩𝑩}^{\mathrm{H}}$ je vždy pozitivně definitní:

${𝑩}^{-1}={𝑩}^{\mathrm{H}}({𝑩𝑩}^{\mathrm{H}}{)}^{-1}$

Metoda nejmenších čtverců

Soustavy $𝑨𝒙=𝒃$ s pozitivně definitní maticí soustavy se v aplikacích objevují poměrně často. Například normálové rovnice v lineárních úlohách nejmenších čtverců mají tento tvar a ostatně i vedly k objevu Choleského rozkladu.

Může se také stát, že matice $𝑨$ pochází z energetického funkcionálu, který musí být z fyzikálních důvodů kladný. Podobný případ často nastává při numerickém řešení parciálních diferenciálních rovnic .

Nelineární optimalizace

Nelineární vícerozměrné funkce mohou být minimalizovány přes jejich parametry pomocí variant Newtonovy metody nazývané kvazi-Newtonovy metody. Při $k$-té iteraci se postupuje k řešení ve směru ${𝒑}_k$ definovaným řešením soustavy ${𝑩}_k{𝒑}_k=-{𝒈}_k$, kde ${𝒈}_k$ je gradient a ${𝑩}_k$ je pozitivně definitní aproximace Hessovy matice.

Další aplikace

Mimo matematiku se Choleského rozklad využívá také v ekonometrickém výzkumu makroekonomických vztahů. V tzv. vektorových autoregresních modelech (VAR) se určuje pořadí, ve kterém se endogenní proměnné navzájem ovlivňují.

Kromě toho se také používá v metodě Monte Carlo k přenesení předem určených korelací do nezávisle generovaných sekvencí náhodných čísel jako diskretizace náhodných procesů.

Implementace

V jazyku C lze výpočet rozkladu zapsat následovně:

for (c=0; c<n; c++) {
  for (sum=0, i=c-1; i>=0; i--)
    sum += sqr(L[c][i]);
  L[c][c] = sqrt(A[c][c] - sum);
  for (r=c+1; r<n; r++) {
    for (sum=0, i=c-1; i>=0; i--)
      sum += L[r][i]*L[c][i];
    L[r][c] = (A[r][c] - sum) / L[c][c];
  }
}

Implementace v programovacích knihovnách

• Programovací jazyk C : Vědecká knihovna GNU poskytuje několik implementací Choleského rozkladu.
• Systém počítačové algebry Maxima : funkce cholesky počítá Choleského rozklad.
• Systém numerických výpočtů GNU Octave poskytuje několik funkcí pro výpočet, aktualizaci a aplikaci Choleského rozkladu.
• Knihovna LAPACK poskytuje vysoce výkonnou implementaci Choleského rozkladu, která je přístupná z Fortranu, C a většiny jazyků.
• V Pythonu provádí funkce cholesky z modulu numpy.linalg Choleského rozklad.
• V Matlabu dává funkce chol Choleského rozklad. Všimněte si, že chol standardně vrací Choleského faktor $𝑼$ v horním trojúhelníkovém tvaru, tj. počítá rozklad $𝑨={𝑼}^{\mathrm{H}}𝑼$. Lze předat příznak, aby se místo toho použil dolní trojúhelníkový faktor.
• V R dává funkce chol Choleského rozklad.
• V Julia poskytuje funkce cholesky ze standardní knihovny LinearAlgebra Choleského rozklad.
• V Mathematice lze na matici aplikovat funkci „CholeskyDecomposition“.
• V C++ podporuje tento rozklad několik knihoven lineární algebry:
  • Armadillo (knihovna C++) poskytuje příkaz chol k provedení Choleského rozkladu.
  • Knihovna Eigen poskytuje Choleského rozklady pro řídké i husté matice.
  • V balíčku ROOT je k dispozici třída TDecompChol .

• V Analytica poskytuje funkce Decompose Choleského rozklad.
• Knihovna Apache Commons Math má implementaci Šablona:Wayback, kterou lze použít v Javě, Scale a jakémkoli jiném jazyce JVM.

Odkazy

Literatura

• Šablona:Citace monografie
• Šablona:Citace monografie
• Šablona:Citace monografie
• Šablona:Citace elektronické monografie
• Šablona:Citace elektronické monografie

Související články

• Definitnost
• Hermitovská matice