Pythonプログラミング(ステップ7・配列・線形代数・行列の分解)
このページでは、NumPyライブラリの線形演算機能を用いる例として、各種の行列の分解についてまとめておく。
行列を行列の積に分解する代表的な操作
LU分解とコレスキー分解
連立一次方程式 $$ A \boldsymbol{x} = \boldsymbol{b} $$ を解く際に、もしも行列$A$が、下三角行列$L$と上三角行列$U$の積 $$ A = L U $$ に分解できるとすると(これはLU分解と呼ばれる)、解くべき方程式は $$ L U \boldsymbol{x} = \boldsymbol{b} $$ となるから、 $$ \boldsymbol{y} = U \boldsymbol{x} $$ と未知変数を置き直せば、まず、 $$ L \boldsymbol{y} = \boldsymbol{b} $$ の解を求め、次に $$ U \boldsymbol{x} = \boldsymbol{y} $$ を解けば、最終的に解 $\boldsymbol{x}$ を求めることができる。
問題を余計に面倒にしただけのようにも思えるが、三角行列で与えられた連立方程式は、掃き出し法でとても簡単に解くことができるので、 一旦LU分解が出来たら、異なる$\boldsymbol{b}$に対して解を得るのは容易である。
また、三角行列の対角成分はその固有値が並んでいるので、その行列式は対角成分の積に等しいことから、$|A| = | LU | = |L| |U|$ によって、LU分解を経由して$A$の行列式を求めることができる。
$A$が正定値対称行列(固有値が全て正であるような対称行列)ならば、下三角行列とその転置のみを使って $$ A = L L^\top $$ と一意に分解できる。これは、特にコレスキー分解と呼ばれている。
QR分解
任意の正方行列$A$は、直交行列$Q$と、上(右)三角行列$R$の積 $$ A = QR $$ に分解することができる(QR分解)。 $A$の基底と、直交基底が、三角行列によって媒介される格好の式になっている。
行列$A$を列ベクトルの並び $$ A = \left( \boldsymbol{a}_1 \boldsymbol{a}_2 \cdots \boldsymbol{a}_n \right) $$ とみなして、これらのベクトルからグラム・シュミットの方法で、正規直交系を構成する手順を考える。
$$ \begin{eqnarray} \boldsymbol{q}_1 & = & \frac{\boldsymbol{a}_1}{\| \boldsymbol{a}_1 \|} \\ \boldsymbol{v}_2 & = & \boldsymbol{a}_2 - (a_2 \cdot q_1) q_1, \;\; \boldsymbol{q}_2 = \frac{\boldsymbol{v}_2}{\| \boldsymbol{v}_2 \|} \\ \vdots & &\\ \boldsymbol{v}_j &=& \boldsymbol{a}_j - (a_j \cdot q_1) q_1 - (a_j \cdot q_2) q_2 - \cdots - (a_j \cdot q_{j-1}) q_{j-1}, \;\; \boldsymbol{q}_j = \frac{\boldsymbol{v}_j}{\| \boldsymbol{v}_j \|} \\ \end{eqnarray} $$
同じ式を並べ替えて $$ \begin{eqnarray} \boldsymbol{a}_1 & = & \| \boldsymbol{a}_1 \| \boldsymbol{q}_1 \\ \boldsymbol{a}_2 & = & (\boldsymbol{a}_2 \cdot \boldsymbol{q}_1) \boldsymbol{q}_1 + \| \boldsymbol{v}_2 \| \boldsymbol{q}_2 \\ \vdots & & \\ \boldsymbol{a}_j & = & (a_j \cdot q_1) q_1 - (a_j \cdot q_2) q_2 - \cdots - (a_j \cdot q_{j-1}) q_{j-1} + \| \boldsymbol{v}_j \| \boldsymbol{q}_j \\ \end{eqnarray} $$ これを行列の形 $$ \left( \boldsymbol{a}_1 \boldsymbol{a}_2 \cdots \boldsymbol{a}_n \right) = \begin{pmatrix} \boldsymbol{q}_1 & \boldsymbol{q}_2 & \cdots & \boldsymbol{q}_n \\ \end{pmatrix} \begin{pmatrix} \| \boldsymbol{a}_1 \| & \boldsymbol{a}_2 \cdot \boldsymbol{q}_1 & \cdots & \boldsymbol{a}_n \cdot \boldsymbol{q}_1 & \\ 0 & \| \boldsymbol{v}_2 \| & & \boldsymbol{a}_n \cdot \boldsymbol{q}_2 & \\ \vdots \\ 0 & 0 & \cdots & \| \boldsymbol{v}_n \| \\ \end{pmatrix} $$ にまとめてみると、グラムシュミットの直交化法によって、$A$が直交行列と上三角行列の積に分解されていることが確認できる。
対角行列を含む分解(対角化)
データ処理にしばしば登場する相関行列や、電気回路、制御、化学プロセス、物理現象等に関係してしばしば登場する行列の多くは対称行列$S$である。 実対称行列の固有値は実数で、その固有値を対角要素に並べた対角行列$\Lambda$と、直行行列$P$を使って $$ S = P \Lambda P^\top $$ のように分解することができる。 ここで、$P^{-1} = P^\top$である。
実対称行列を複素数に拡張したエルミート行列$U$は $$ U = P \Lambda P^* $$ のように分解できる。 ここで、$P^*$は$P$のエルミート共役(転置を行い、要素の複素共役を取った行列)である。
$n$次の正方行列$A$を構成する$n$個の固有ベクトル全てが互いに独立ならば、 固有ベクトル(列ベクトル)を並べて構成した行列$P$を使って、 $$ A = P \Lambda P^{-1} $$ と分解できる。 ここでも$\Lambda$は、対角要素に固有値を並べた対角行列である。
特異値分解
任意の$m$行$n$列の行列$A$は、$m$行$m$列のユニタリー行列 $U$ と $n$行$n$列のユニタリー行列$V$のエルミート共役$V^*$を使って $$ A = U \Sigma V^* $$ のように分解できる。 ここで、行列$A$から構成される列(あるいは行)ベクトルの独立な成分の数(ランク)を$r$とすれば、 $\Sigma$は、サイズ$r$の対角ブロックと、0要素で構成されるブロックで構成される$m$行$n$列の行列である: $$ \Sigma = \left( \begin{array}{@{}c|c@{}} \begin{matrix} \lambda_1 & 0 & & 0 & \\ 0 & \lambda_2 & & 0 & \\ & & \ddots & & \\ 0 & 0 & & \lambda_r & \\ \end{matrix} & O \\ \hline O & O \\ \end{array} \right) $$
ユニタリー行列(直交行列)$U$と$V$はベクトルの回転や反転、対角行列$\Sigma$は伸長や短縮を表すことを思い出すと、任意の線形変換は、回転・反転の後に、軸に沿った伸縮操作を施し、再び回転・反転を行うという3ステップに分解できる、とも解釈できる。
NumPyの機能
NumPyの線形代数ライブラリ(numpy.linalg)には、行列の分解がいくつか実装されている。コードの冒頭で
import numpy as np
が宣言されていることを前提に、表にまとめておく。
関数 | 機能 |
---|---|
L = np.linalg.cholesky(A) |
対称行列Aをコレスキー分解した行列を返す |
Q,R = np.linalg.qr(A) |
AのQR分解 |
w,P = np.linalg.eigh(U) |
対称行列またはエルミート行列の固有値の一次元配列wと 固有ベクトルを並べた行列(すなわち変換行列)Pを返す。 †固有ベクトルは縦(列)ベクトルとして返される |
w,P = np.linalg.eig(A) |
正方行列の固有値wと対応する固有ベクトルを並べた行列Pを返す |
U,s,V = np.linalg.svd(A,full_matrices=False) |
行列を特異値分解した結果を返す。 sは対角要素の一次元配列。 ここで返されるVは$V^*$。 |
対角成分が一次元配列s
で与えられているとき
np.diag(s)
で(非対角要素の0が補われた)対角行列を得ることができる。