線形代数
QR分解; QR Decomposition

概要

---

QR分解はある行列\(A\)を直交行列\(Q\)と上三角行列\(R\)を用いて\(A=Q R\)として分解する手法である。
QR分解することで固有値が計算しやすくなる(収束速度と安定性が良くなる)という利点がある。

Gram-Schmidt法による\(Q, R\)は以下のように与えられる。
\( \quad A = \left( \boldsymbol{a_1} \boldsymbol{a_2} \cdots \boldsymbol{a_N} \right) , \quad \boldsymbol{u_i} = \boldsymbol{a_i} - \displaystyle \sum_{j=1}^{i-1} \boldsymbol{u_j} \frac{ \boldsymbol{a_i} \cdot \boldsymbol{u_j} }{ |\boldsymbol{u_j}|^2 }, \quad \boldsymbol{e_i} = \frac{\boldsymbol{u_i}}{|\boldsymbol{u_i}|} \)
\( \quad Q = \left( \boldsymbol{e_1} \boldsymbol{e_2} \cdots \boldsymbol{e_N} \right) , \quad R = \begin{pmatrix} \boldsymbol{a_1} \cdot \boldsymbol{e_1} & \boldsymbol{a_2} \cdot \boldsymbol{e_1} & & \cdots & & \boldsymbol{a_N} \cdot \boldsymbol{e_1} \\ 0 & \boldsymbol{a_2} \cdot \boldsymbol{e_2} & & & & \\ \vdots & \ddots & \ddots & & & \vdots \\ & & & \boldsymbol{a_i} \cdot \boldsymbol{e_i} & & \\ \vdots & & & \ddots & \ddots & \\ 0 & \cdots & & \cdots & 0 & \boldsymbol{a_N} \cdot \boldsymbol{e_N} \end{pmatrix} \)

行列クラス本体はこちら

QR分解の性質

---

QR分解で求められた\(Q, R\)は以下のような性質を持つ。
\( \quad A = QR \)
\( \quad R = Q^{T} A \)


ソースコード

---

namespace Algebra {
    /// <summary>行列クラス</summary>
    public partial class Matrix {
        /// <summary>QR分解</summary>
        public void QRDecomposition(out Matrix orthogonal_matrix, out Matrix triangular_matrix) {
            if(!IsSquare(this)) {
                throw new InvalidOperationException();
            }

            orthogonal_matrix = new Matrix(Size, Size);
            triangular_matrix = new Matrix(Size, Size);

            int i, j, n = Size;

            Vector[] e = new Vector[n], u = new Vector[n];
            for(i = 0; i < Size; i++) {
                e[i] = Vector.Zero(n);
            }

            for(i = 0; i < n; i++) {
                Vector ai = Vertical(i);

                u[i] = ai.Copy();

                for(j = 0; j < i; j++) {
                    u[i] -= u[j] * Vector.InnerProduct(ai, u[j]) / u[j].SquareNorm;
                }

                e[i] = u[i].Normal;

                for(j = 0; j <= i; j++) {
                    triangular_matrix[j, i] = Vector.InnerProduct(ai, e[j]);
                }

                for(j = 0; j < n; j++) {
                    orthogonal_matrix[j, i] = e[i][j];
                }
            }
        }
    }
}

関連項目

---

行列
ガウスの消去法
LU分解
行列式
トレース
固有値・固有ベクトル