JA EN
体系線形代数
·★ 会員·17分で読めます

行列分解ツアー — LU・QR・コレスキー・SVDの使い分け

行列分解は「同じ変形を、扱いやすい部品の掛け算に書き直す」操作です。解く(LU)・半額で解く(コレスキー)・直交化して最小二乗(QR)・壊れた行列でも扱う(SVD)の4つを、なぜ逆行列を作ってはいけないのかから始めて、選び方の一枚表とnumpyの実装まで一気につなぎます。

対象textタスクmath

逆行列を作らずに「解く」

学校で Ax=bAx = b を解けと言われたら、まず x=A1bx = A^{-1}b と書きます。ところが数値計算の現場では、逆行列を作る関数はほとんど使われません。ライブラリのドキュメントにも「inv の代わりに solve を使え」と書いてあります。

理由は2つあります。1つは、逆行列を全部求めるのは、方程式を1回解くよりはっきり高くつくこと。もう1つは、逆行列を経由すると誤差が余計に混ざることです。求めたいのは xx という1本のベクトルなのに、わざわざ n×nn \times n の表を丸ごと作ってから掛け算する。ゴールまで直線で行けるのに、いったん遠回りしているようなものです。

では何をするのか。行列分解です。整数の 60 を 22352^2 \cdot 3 \cdot 5 と書き直せば約数がひと目で分かるように、行列も「扱いやすい部品の積」に書き直すと、そこから先が急に簡単になります。この記事で扱うLU・コレスキー・QR・SVDは、どれも同じ発想の別バージョンです。違うのは「元の行列にどんな性質を仮定するか」と「その見返りにどれだけ安く・頑丈に済むか」だけです。

前作LoRAとRAGを支える線形代数では、行列を「空間の変形装置」として見ました。今回はその変形を分解して運ぶ話です。

扱いやすい部品は、実は2種類しかない

分解の行き先は、驚くほど少ない種類しかありません。三角行列と、直交行列です。

三角行列は、対角より下(または上)が全部ゼロの行列です。なぜこれが嬉しいかというと、方程式がドミノ倒しで解けるからです。

Ly=byi=1ii(bij<iijyj)L\boldsymbol{y} = \boldsymbol{b} \quad\Longrightarrow\quad y_i = \frac{1}{\ell_{ii}}\Big(b_i - \sum_{j<i} \ell_{ij}\, y_j\Big)
(1)

この式が言っているのは、ii 番目の答えは、それより前に確定した答えを代入して引き算するだけで決まるということです。1行目は変数が1つしかないので即決まる。2行目は1行目の答えを代入すれば1つしか残らない。以下同様。連立方程式なのに、実質は上から順の代入作業です。この手順を前進代入(下三角)・後退代入(上三角)と呼び、コストは n2n^2 に比例します。n3n^3 ではありません。

直交行列は、列同士が直角で、どの列も長さ1の行列です。こちらの嬉しさは2つあります。

QQ=I,Qx=xQ^\top Q = I,\qquad \|Q\boldsymbol{x}\| = \|\boldsymbol{x}\|
(2)

左は「転置がそのまま逆行列になる」、右は「掛けても長さが変わらない」と言っています。逆行列を求める計算が要らないのが1つ目の嬉しさ。長さを変えないので、誤差を掛け算のたびに膨らませないのが2つ目です。数値計算で直交行列が偏愛されるのは、この「誤差を増やさない」性質のためです。

直交とは要するに内積がゼロということです。下の図で2本のベクトルを回して、直角のとき内積がちょうど0になるのを確かめてください。分解で出てくる QQ は、この「互いに直角で長さ1」を全列で満たした行列です。

FIG 12本のベクトルを回すと内積が正→0→負と変わる。ちょうど直角のときだけ0。直交行列とは、この関係をすべての列の組で満たす行列のこと

つまり行列分解とは、元の行列を「ドミノ倒しで解ける部品」と「誤差を増やさない部品」の積に書き直す作業です。以下の4つの分解は、この2種類をどう組み合わせるかのバリエーションにすぎません。

LU分解 — 消去法を記録に残す

最初の1つは、中学高校でやった加減法(掃き出し法)そのものです。「1行目を2倍して2行目から引く」という操作を繰り返して、下三角部分をゼロにしていく。あれをコンピュータにやらせて、消去に使った倍率を捨てずに保存したものがLU分解です。

PA=LUPA = LU
(3)

LL は対角が1の下三角行列(消去に使った倍率が入っている)、UU は消去後に残った上三角行列、PP は行を入れ替えた記録です。この式は、「行を適当に並べ替えれば、どんな正方行列も下三角×上三角に書ける」と言っています。

PP が要るのは、割り算の分母(ピボット)に0や極端に小さい数が来ると破綻するからです。毎ステップで列の中の絶対値最大の行を先頭に持ってくる——これを部分ピボット選択と呼び、実用のLU分解には必ず付いています。教科書の消去法が実装で壊れる原因のほとんどはここで、その対策が PP 一文字に凝縮されています。

分解さえ済めば、解くのは2段階です。LUx=PbLUx = Pb を、まず Ly=PbLy = Pb を前進代入で、次に Ux=yUx = y を後退代入で解く。分解が n3n^3 に比例するのに対し、この2段階は n2n^2 です。

ここに実務上の最大の利点があります。同じ AA に対して bb だけが何度も変わる場合、分解は1回で済むのです。100本の右辺があっても、重い計算は最初の1回だけ。逆行列を作って100回掛けるより速く、しかも誤差が小さい。「inv ではなく solve」の中身は、要するにこれです。

コスト感も押さえておきます。LU分解に必要な乗除算はおよそ 23n3\tfrac{2}{3}n^3 回です。nn が2倍になると8倍、10倍になると1000倍。n=1000n=1000 なら数億回で、いまのCPUなら1秒かかりません。n=10,000n=10{,}000 ではその千倍になり、ここで初めて「密行列のまま解く」という設計そのものを見直す段階に入ります。一方の前進・後退代入は n2n^2 回で、n=1000n=1000 でも100万回どまり。重いのは分解だけで、解くのはほとんどタダという非対称が、この2つの指数の差です。分解を使い回す価値は、ここから出てきます。

LUにはもう1つ副産物があります。行列式です。detA=±iuii\det A = \pm\prod_i u_{ii}——つまりUU の対角に並んだ数を全部掛け、行を入れ替えた回数の偶奇で符号を決めるだけ、という式です。定義どおりの余因子展開は項数が n!n! に比例し、n=20n=20 でも現実的な時間では終わりませんが、LUを経由すれば分解のついでに手に入ります。ライブラリの行列式関数の中身は、たいていこれです。

行列に良い性質があるなら、その分だけ安く済ませたい。最も見返りが大きい性質が対称正定値です。

この先にあるもの

§

ここから先は会員限定です

解説記事371本・教科書26章・学生モード48単元・論文精読6本が、月額¥490ですべて読み放題になります。新しい解説は毎日3本ずつ増えます。いつでも解約でき、解約後も期間の終わりまで読めます。

会員の方はログインすると続きが表示されます

コメント

コメントにはログインが必要です