QR分解

QR分解（キューアールぶんかい、）とは、m × n 実行列 Aを、 m 次直交行列 Q と m × n 上三角行列 R との積への分解により表すこと、またはそう表した表現をいう。このような分解は常に存在する。

QR分解は線型最小二乗問題を解くために使用される。また、固有値問題の数値解法の1つであるQR法の基礎となっている。

定義

正方行列
すべての実正方行列 Aは直交行列Qと上三角行列(別名右三角行列)Rを用いて

と分解できる。もしAが正則ならば、Rの対角成分が正になるような因数分解は一意に定まる。

もしAが複素正方行列ならば、Qがユニタリ行列 (つまり)となるような分解A = QRが存在する。

もしAがn個の線形独立な列を持つなら、Qの最初のn列はAの列空間の正規直交基底をなす。より一般的に、1 ≤ k ≤ nの任意のkについて、Qの最初のk列はAの最初のk列の線型包をなすL. N. Trefethen and D. Bau, Numerical Linear Algebra (SIAM, 1997).。Aの任意の列kがQの最初のk列にのみ依存するということは、Rが三角行列であることから明らかである。

矩形行列
より一般的に、m ≥ nである複素m×n行列Aを、m×mユニタリ行列Qとm×n上三角行列Rに分解することができる。m×n上三角行列の下から(m−n)行はすべてゼロであるため、Rや、RとQ両方の分割を簡単に行うことができる。

ここで、R1はn×n上三角行列、0は(m − n)×n零行列、Q1はm×n行列、Q2はm×(m − n)行列で、Q1とQ2は両方直交する列を持つ。

Q1R1をはAの薄い(thin)QR分解と呼び、 Trefethen & Bauは軽減(reduced)QR分解と呼んでいる。もしAが最大階数nであり、R1の対角成分を正にするならば、R1とQ1は一意に定まる。しかし一般的にQ2はそうではない。R1はA* A (Aが実行列の場合ATAに等しい)のコレスキー分解の上三角部分に等しい。

QL・RQ・LQ分解
同様に、Lを下(lower)三角行列として、QL、RQ、LQ分解を定義することができる。

QR分解の計算
QR分解を計算する手法として、グラム・シュミット分解、ハウスホルダー変換、ギブンス回転などがある。それぞれ利点と欠点がある。

グラム・シュミットの正規直交化法の使用

グラム・シュミットの正規直交化法を最大階数行列の列に適用することを考える。内積 (複素ベクトルの場合 )とする。

射影の定義より、

したがって、

ここでを新しく計算された正規直交基底上に表すことができ、であるから、

これは行列の形に書くことができ、

ただし、

例

の分解を考える。

正規直交行列に対して、


が常に成り立つから、グラム・シュミット法により以下の手順でを計算できる。


したがって、


RQ分解との関係
RQ分解は行列Aを上三角行列R(別名右三角)と直交行列Qに変換する。QR分解との違いはこれらの行列の順番だけである。

QR分解はグラム・シュミットの正規直交化法をAの最初の列から最後の列の順に適用する。

RQ分解はグラム・シュミットの正規直交化法をAの最後の行から最初の行の順に適用する。

利点と欠点

グラム・シュミットの正規直交化法は本来、数値的に不安定である。射影の応用として直交化との幾何学的な類似性があるが、直交化自体は数値的な差異が生じやすい。しかしながら、実装が簡単という大きな利点があり、外部線形代数ライブラリが利用できない場合のプロトタイピングに便利なアルゴリズムである。

ハウスホルダー変換の使用

ハウスホルダー変換 (またはハウスホルダー鏡映)はベクトルを取り、平面または超平面に関する鏡映をする変換である。この演算はm×n行列 (m ≥ n)のQR変換の計算に使うことができる。

Qは一つの座標を除いたすべての座標が未知でもベクトルを鏡映するために使用できる。

スカラαに対してを満たすようなの任意の実m次元列ベクトルを考える。もしこのアルゴリズムが浮動小数点演算を用いて実装されている場合、桁落ち（による誤差の拡大）を防ぐため、行列Aの最終的な上三角部分においてすべての要素が0である後のピボット座標に対し、αをのk番目の座標の逆符号とする。複素行列の場合には、

として、さらに以下のQの導出において転置を共役転置に読み替えること。

ここで、をベクトル(1 0 … 0)T、||·||をユークリッド距離、をm×m単位行列とし、


あるいは、が複素行列ならば、

とおく。

はm×mハウスホルダー行列であり、


これによりm×n行列Aを上三角の形に漸次変換できる。まず、xの最初の列を選んで得られるハウスホルダー行列Q1にAを乗算する。この結果、行列Q1Aは左の列が(最初の行を除き)ゼロになる。


この操作をA′ (Q1Aから最初の列と最初の行を除いたもの)に繰り返すと、ハウスホルダー行列Q′2が得られる。Q′2はQ1より小さいということに注意すること。A′の代わりにQ1Aで計算したいため、A′を左上に拡張し、ひとつの1を埋める必要がある。つまり、一般的には

とする。
回このプロセスを繰り返すと、のとき、

は上三角行列である。そこで、

とすると、
はのQR分解である。

この鏡映変換を用いた計算方法は先述のグラム・シュミット法よりも数値的安定性がある。

下表にサイズnの正方行列を仮定したときの、ハウスホルダー変換によるQR分解のk番目のステップにおける計算量を示す。
 演算 k番目のステップにおける計算量 乗算  加算  除算  平方根 

これらの数をn − 1ステップ(サイズnの正方行列であるため)まで合計して、このアルゴリズムの(浮動小数点演算の観点からの)複雑さは

と表せる。

例

の分解を考える。

まず、行列Aの最初の列、ベクトルを、に変換する鏡映を見つける必要がある。 

今、

とする。

ここで、
 であり、

であるから、
  and であり、

である。

今、

を見てみると、すでにほぼ三角行列である。あとは(3, 2)要素を零にするだけでよい。

(1, 1)における小行列を取り、同じプロセスを

に再び適用する。

先述のメソッドと同様にして、このプロセスの次のステップが正しく動作するために、1で直和を取ることにより、ハウスホルダー変換

を得る。

今、

または、有効数字四桁で、

を得た。

行列Qは直交行列であり、Rは上三角行列であるため、A = QRは求めるべきQR分解である。

利点と欠点

ハウスホルダー変換の使用は、R行列のゼロを生成するメカニズムに鏡映を利用しているため、最もシンプルでかつ数値的に安定したQR分解アルゴリズムである。しかしながら、新しい零要素を生成する毎回の鏡映変化において行列QとR両方の行列全体を書き換えるため、ハウスホルダー変換は必要とするメモリ帯域幅が多く、並列化できない。（ただし鏡映変換のブロック化に基づくブロック化ハウスホルダ変換を使えばメモリ帯域幅への要求を低減することができる）

ギブンス回転の使用
QR分解はギブンス回転を使用しても計算できる。各回転により行列の亜対角要素がゼロになり、R行列を構成できる。すべてのギブンス回転を結合することで直交行列Qを構成できる。

実際には、行列全体を構成して乗算をするようなギブンス回転は行われない。代わりに、疎な要素を計算するような無駄な計算をしない、疎なギブンス行列乗算と同等なあるギブンス回転の手順が採られる。そのギブンス回転の手順は少しの非対角成分をゼロにするだけで済み、ハウスホルダー変換よりも容易に並列化できる。

例

の分解を考える。

まず、左下隅の要素、をゼロにする回転行列を構成する必要がある。この行列はギブンス回転で求めることができる。まずベクトルを、X軸に沿って回転させる。このベクトルは角度を持つ。直交ギブンス回転行列を次のように作る。

ここでの結果は要素がゼロである。

同様にしてそれぞれ非対角要素・要素がゼロであるようなギブンス行列・を構成し、三角行列を構成する。直交行列はすべてのギブンス行列の積で表される。したがって、であり、QR分解はである。

利点と欠点

ギブンス回転によるQR分解は、アルゴリズムを完全に動かすのに必要な行の順序を決定するのが簡単ではないため、実装に最も手間がかかる。しかしながら、新しいゼロ要素がゼロになる予定の要素の行(i)とその上の行(j)にしか影響しないという特筆すべき利点がある。これにより、ギブンス回転アルゴリズムはハウスホルダー変換手法よりも帯域幅効率が良く、容易に並列化できる。

行列式や固有値の積との関係
QR分解を正方行列の行列式の絶対値を求めるのに利用できる。ある行列がと分解できるとする。このとき

である。

Qはユニタリであるため、である。したがって、をRの対角要素とすると、

となる。

さらに、行列式は固有値の積に等しいため、をの固有値とすると、


となる。

QR分解の定義を非正方行列に導入し、固有値を特異値に置き換えることで、上記性質を非正方行列に拡張することができる。

非正方行列AのQR分解を


とする。ただし、は零行列、はユニタリ行列。

特異値分解と行列式の性質から、をの特異値として、

との特異値は同じであるが、複素固有値が異なる場合があることに注意すること。しかしながら、Aが正方ならば、下記は真である。

結論として、QR分解を使うことによって行列の固有値や特異値の積を効率よく計算することができる。

列のピボット
ピボットQRは列のピボットの新しいステップにおいて、それぞれ初めに残りの列で最も大きいものを取るという点で、通常のグラム・シュミット法とは異なっている。したがって、置換行列Pを次のように導入する。

列のピボットはAが(ほぼ)階数落ちである、またはその疑いがある場合に便利である。また、数値的精度を向上させることもできる。通常、Rの対角成分が非増加、つまりとなるようにPを選ぶ。この手段により特異値分解よりも低い計算コストでAの(数値的な)階数を求めることができ、いわゆるの基礎となっている。

線形逆問題への利用
行列の逆行列を直接求めるのに比べ、QR分解を利用した逆問題の解法は、条件数が減少していることからも分かるように、数値的に安定しているParker, Geophysical Inverse Theory, Ch1.13.。

次元がで階数がであるような行列に対して、劣決定()線形問題を解くためには、まずの転置行列のQR分解を求める。ただし、Qは直交行列(つまり、)であり、Rはという特殊な形である。ここで、は正方右三角行列、零行列は次元である。計算すると、この逆問題の解を次のように表すことができる。

ここで、はガウスの消去法で計算でき、はを用いることで直接計算できる。後者の手法の方が数値的精度が高く、計算量も少ないという利点がある。

ノルムを最小にするような過決定()問題の解を求めるためには、まずのQR分解を求める。を直交行列全体のうち最初の列を含む行列、を先述の通りに置くと、この問題の解はと表せる。劣決定の場合と同様に、の逆行列を直接計算しなくても、を用いることで早く正確にを求めることができる。(とは数値ライブラリによっては高速な(economic)QR分解として実装されている。)

一般化
はQR分解を半単純リー群に一般化している。

 脚注 

 参考文献 
和文





英文






 .

 関連項目 
 QR法
 極分解
 固有値分解
 LU分解
 特異値分解

外部リンク
Online Matrix Calculator 行列のQR分解の実行
LAPACK users manual QR分解の計算のサブルーチンの詳細
Mathematica users manual QR分解の計算のルーチンの詳細と例
ALGLIB LAPACKのC++、C#、Delphiなどへの移植
Eigen::QR QR分解のC++での実装

Category:数値線形代数
Category:行列
Category:行列の分解
Category:数学に関する記事