🗂

QR分解とハウスホルダー変換

2021/03/21に公開

(個人の勉強ノートです。正確さは保証しません。あしからず。)

An 次元の正則な正方行列とする。

A=QR

ここでQ は直交行列、R は上三角行列。

※ より一般にはQR分解は m\times n\ (m\ge n) の行列に対して定義できるが、簡単のため正方行列に限定して勉強。

シュミットの直交化による方法

シュミットの直交化とは、n本の線形独立なベクトル a_1,\cdots,a_n から正規直交基底 u_1,\cdots,u_n を構成する手法である。
つまり行列 A=(a_1\cdots a_n) から直交行列 Q=(u_1\cdots u_n) を作る手法といえる。
これを使って、A=QR の形で表せる。これがQR分解で、Rは上三角行列となる。

シュミットの直交化:

\begin{array}{ll} u_1&=a_1/\|a_1\| \\ u_2&=\frac{u_2'}{\|u_2'\|},\quad u_2'=(I-u_1u_1^T)a_2 \\ u_3&=\frac{u_3'}{\|u_3'\|},\quad u_3'=(I-u_1u_2^T-u_2u_2^T)a_3\\ \cdots \end{array}

以上から分かるように、

{\rm span}\{u_1,\cdots,u_k\}={\rm span}\{a_1,\cdots,a_k\}\quad(1\le k\le n)

である。よって、a_k\{u_1,\cdots,u_k\} の線形和で表現できることになる。すなわち、

\begin{array}{ll} a_1&=r_{11}u_1 \\ a_2&=r_{12}u_1 + r_{22}u_2 \\ & \vdots \\ a_k&=r_{1k}u_1 + r_{2k}u_2 + \cdots + r_{kk}u_k \\ & \vdots \end{array}

と表現できる。これを行列として書き表すと次式となる(QR分解)。

A=(u_1\cdots u_n)\left(\begin{array}{cccc} r_{11} & r_{12} & \cdots & r_{1n} \\ & r_{22} & \cdots & r_{2n} \\ & & \ddots & \vdots \\ & & & r_{nn} \\ \end{array}\right) = QR

具体的な r_{ij} を求めるには、

a_j=r_{1j}u_1 + r_{2j}u_2 + \cdots + r_{jj}u_j=\sum_{i=1}^j r_{ij}u_j

の両辺に u_i との内積を取れば良い。

(u_i, a_j)=\sum_{i=1}^j r_{ij}(u_i, u_j)= r_{ij}\quad(\because (u_i,u_j)=\delta_{ij})

ハウスホルダー変換による方法

シュミットの直交化による方法は簡便であるが、数値的に不安定である。従って、実用的には数値的により安定なハウスホルダー変換による方法が用いられる。

ハウスホルダー変換とは

法線ベクトル n の(原点を通る)超平面に対して対称に点を写す鏡映変換のこと。この変換行列を H で表す。(※ \|n\|=1

Hx = x - 2(n, x)n = (I-2nn^T)x

これは直交行列であり、かつ対称行列である。

対称性:

H^T=(I-2nn^T)^T=(I^T-2nn^T)=H

直交性:

\begin{array}{ll} H^TH&=H^2\\ &=I^2-4nn^T+4(nn^T)(nn^T) \\ &=I-4nn^T+4n(n^Tn)n^T \\ &= I \end{array}

単なる鏡映変換に「ハウスホルダー変換」はちょっと大仰な名前に感じられるが(※ 個人の感想です)、上記のような特筆すべき性質が大仰な名前に値するということかもしれない。実際、これらの性質が以下で大活躍している。簡単なのに名前でドヤれるので勉強の投資対効果が高い(数学をそういうことにつかってはいけない)。

ハウスホルダー変換によるQR分解の骨子

A_0=A とおいて

A_k=H_kA_{k-1}

という漸化式で A に次々とハウスホルダー変換を施すと、A_n は次式となる。

\begin{array}{ll} A_n&=H_nA_{n-1} \\ &=H_nH_{n-1}A_{n-2} \\ &=\cdots \\ &=H_nH_{n-1}\cdots H_1 A \end{array}

この A_n が上三角行列 R となるようにハウスホルダー変換の列 H_k を定めることができたとする(この定め方は後述)。すると次式が得られる。

H_nH_{n-1}\cdots H_1 A=R\\

これを A について解く。ここでハウスホルダー変換が対称行列かつ直交行列であるという性質が効いてくる。

\begin{array}{lll} A&=(H_nH_{n-1}\cdots H_1)^{-1}R \\ &=H_1^{-1}H_2^{-1}\cdots H_n^{-1} R \\ &=H_1^TH_2^T\cdots H_n^T R &\because H_kは直交行列 \\ &=H_1H_2\cdots H_n R &\because H_kは対称行列 \\ &=QR \end{array}

直交行列の積は直交行列であるから、Q は直交行列である。つまりこれはQR分解になっている。

続いて、ハウスホルダー変換の列 \{H_k\} の定め方について述べるが、その前に幾つか準備がある。

準備1: xy に写すハウスホルダー変換

まず準備として、2点 x, y が与えられたときに Hx=y となる H を求めておく。

任意の x, y に対して H が定義できるわけではない。H は直交行列であるから、変換の前後でノルムは保存される。すなわち \|x\|=\|y\| でなくてはならない。

鏡像変換を作る超平面の法線ベクトルは次式で与えられる。

n=\frac{x-y}{\|x-y\|}

よってハウスホルダー行列は次式。

H=I-2\frac{(x-y)(x-y)^T}{\|x-y\|^2}

特にQR分解の文脈では、次のようなハウスホルダー変換を行う。

x=\left(\begin{array}{c}x_1\\x_2\\ \vdots\\ x_n\end{array}\right),\quad y=Hx=\left(\begin{array}{c}\|x\| \\ 0\\ \vdots\\ 0\end{array}\right)

準備2: 部分空間に対するハウスホルダー変換

法線ベクトルが

n=\left(\begin{array}{c}0 \\ \tilde{n}\end{array}\right)

という形をしている超平面によるハウスホルダー変換を考える。(0\in R^k,\ \tilde{n}\in R^{n-k}

\begin{array}{ll} H&=I-2nn^T \\ &= I-2\left(\begin{array}{c}0\\ \tilde{n}\end{array}\right)\left(\begin{array}{cc}0^T&\tilde{n}^T\end{array}\right) \\ &= I-2\left(\begin{array}{cc} 0 & 0 \\ 0 & \tilde{n}\tilde{n}^T \end{array}\right) \\ &= \left(\begin{array}{cc} I & 0 \\ 0 & I-2\tilde{n}\tilde{n}^T \\ \end{array}\right) \\ &= \left(\begin{array}{cc} I & 0 \\ 0 & \tilde{H} \\ \end{array}\right) \end{array}

\{H_k\} の定め方

準備が整ったので、QR変換を導くハウスホルダー変換の列を導出していく。

変換の漸化式を再掲する。

A_k=H_kA_{k-1}

この A_k は次の形となる(そのように H_k を構成する)。

A_k=\left(\begin{array}{cc} R_k & B_k \\ O & C_k \end{array}\right)

ここで R_kk次元の上三角行列、B_kk\times (n-k)行列、C_kn-k次元の正方行列。(k=0A_0=B_0、また k=nA_n=R_n となり上三角行列となることに注目。)

このとき、A_k にハウスホルダー変換 H_{k+1} を作用させて同様の形の A_{k+1} を得たい。

A_{k+1}=H_{k+1}A_k=\left(\begin{array}{cc} R_{k+1} & B_{k+1} \\ O & C_{k+1} \end{array}\right)

同様に R_{k+1}k+1次元の上三角行列、C_{k+1}n-k-1次元の正方行列。

このためには、行列 C_k を次のように列ベクトルが並んでいると考え、

C_k=(c^k_1\ c^k_2\ \cdots\ c^k_{n-k})

先頭の列ベクトル c^k_1 を次のように変換するハウスホルダー変換 \tilde{H}_{k+1} を作れば良い。

\tilde{H}_{k+1}c^k_1=\left(\begin{array}{c}\|c^k_1\| \\ 0\\ \vdots\\ 0\end{array}\right),\quad H_{k+1}=\left(\begin{array}{cc} I & 0 \\ 0 & \tilde{H}_{k+1} \\ \end{array}\right)

これを A_k に作用させると次式となる。

\begin{array}{ll} A_{k+1}&=H_{k+1}A_k \\ &= \left(\begin{array}{cc} I & 0 \\ 0 & \tilde{H}_{k+1} \\ \end{array}\right)\left(\begin{array}{cc} R_k & B_k \\ O & C_k \end{array}\right) \\ &=\left(\begin{array}{cc} R_k & B_k \\ O & \tilde{H}_{k+1}C_k \end{array}\right) \end{array}

以上より、こういったハウスホルダー変換を帰納的に繰り返せばQR変換が得られることが分かる。

参考にした記事

チョットワカル【QR分解】
https://qiita.com/Sharkkii/items/cebba1bc538fdaaa9fb5

Discussion