前書き

記事の概要

今回は行列の低ランク近似など様々な場所で活躍する 特異値分解 について解説していきます。
この記事では

数式による特異値分解のアルゴリズムの解説
プログラム実装例

が含まれています。それぞれ自分に必要な部分を参照してもらえればいいなと思っています。

前提知識

大学1年生で勉強する「線形代数学」の基本的な事項は既知として書いていきます。
キーワードは 直交行列 / 固有値分解 / グラム・シュミットの直交化法 です。
もし記事中でわからない言葉が出てきた場合には、線形代数学の教科書を読んだり他のサイトを見たりして調べてください。

導入

特異値分解って何？

$m \times n$ 行列 $A$ に対して、 $A = U \Sigma V^T$ を満たす $m$ 次直交行列 $U$ 、 $n$ 次直交行列 $V$ 、 $m \times n$ 対角行列 $\Sigma$ を求めることを 特異値分解 と言います。¹
特異値分解はどんな行列に対しても実行できることが知られています。

特異値分解できると何が嬉しいの？

特異値分解によって、その行列に含まれる重要な成分を抽出することができます。重要な成分を見つけることができれば、そうでない成分を削除することで情報量を減らすことができます。

特異値分解の応用先はとても幅広いです。ここでは 行列の低ランク近似 について簡単に紹介します。

行列の低ランク近似 とは、与えられた行列 $A$ をそれよりもランクの低い行列によって近似するという問題です。大きなサイズの行列を「重要な成分を残しながら」小さなサイズの表現形式で近似することで、情報をほとんど損なうことなくデータを圧縮することができます。行列の低ランク近似の記事を書いたので、良ければご覧ください。

理論

基礎的な話

特異値・特異ベクトル

$m \times n$ 行列 $A$ に対して

\begin{eqnarray*}
A {\bf v} &=& \sigma {\bf u} \\
A^T {\bf u} &=& \sigma {\bf v} \\
\end{eqnarray*}

を満たすような $\sigma, {\bf u}, {\bf v}$ をそれぞれ $A$ の 特異値 / 左特異ベクトル / 右特異ベクトル と言います。ここで ${\bf u}$ は $n$ 次元ベクトル、 ${\bf v}$ は $m$ 次元ベクトルです。定義から明らかに、 $A$ が対称行列ならば $A$ の特異値は固有値であり、 $A$ の左・右特異ベクトルは固有ベクトルに一致します。

特異値・特異ベクトルと固有値・固有ベクトルの関係を見るために、特異値に関する定義式

$$
A {\bf v} = \sigma {\bf u} \
$$

の両辺に左から $A^T$ をかけてみます。このとき

$$
A^T A {\bf v} = \sigma A^T {\bf u} = \sigma^2 {\bf v}
$$

が成立します。同様に特異値に関するもう1つの定義式

$$
A^T {\bf u} = \sigma {\bf v}
$$

の両辺に左から $A$ をかけてみます。このとき

$$
A A^T {\bf u} = \sigma A {\bf v} = \sigma^2 {\bf u}
$$

が成立します。この観察から

$A$ の特異値 $\sigma$ は $A A^T, A^T A$ の固有値の(正の)平方根である
$A$ の左特異ベクトル ${\bf u}$ は $A A^T$ の固有ベクトルである
$A$ の右特異ベクトル ${\bf v}$ は $A^T A$ の固有ベクトルである

という性質があることがわかります。

特異値分解のアルゴリズム

$m \times n$ 行列 $A$ に対して特異値分解することを考えます。以下の説明では $m \ge n$ としておきます。( $m \le n$ のときには $A^T$ に対して特異値分解をすれば同じ議論ができるため、一般性を失いません。)

まず、 $A^T A$ に対して固有値分解を実行します。 $A^T A$ は半正定値対称行列なので、 $n$ 次直交行列 $V$ によって対角化できて

$$
V^T A^T A V = \text{diag}(\sigma_1^2, \dots, \sigma_n^2)
$$

が成立します。² ただし、 $\text{diag}(\sigma_1^2, \dots, \sigma_n^2)$ は対角成分を $\sigma_1^2, \dots, \sigma_n^2$ とする $n$ 次対角行列を表しており、 $\sigma_1 \ge \sigma_2 \ge \dots \sigma_k > \sigma_{k+1} = \dots = \sigma_n = 0$ が成立していると仮定します。( $V$ の各列を入れ替えることでこれを満たすようにできます。)

説明の都合上、 $B \equiv A V$ と書くことにします。このとき、 $B^T B$ が対角行列になることから、 $B = ({\bf b_1}, \dots, {\bf b_n})$ の各 ${\bf b_i}$ は互いに直交していることがわかります。式で書けば

$$
\langle {\bf b_i}, {\bf b_j} \rangle = \delta_{ij} \sigma_i \sigma_j
$$

となっています。 $i \le k, j \le k$ なる $i, j$ に限れば $\sigma_i > 0, \sigma_j > 0$ なので

$$
\langle \frac{1}{\sigma_i}{\bf b_i}, \frac{1}{\sigma_j}{\bf b_j} \rangle = \delta_{ij}
$$

と書くことができます。 ${\bf u_i} \equiv \frac{1}{\sigma_i} {\bf b_i} \ (i = 1, \dots, k)$ とすれば

$$
\langle {\bf u_i}, {\bf u_j} \rangle = \delta_{ij}
$$

です。この式は「 $k$ 本のベクトル ${\bf u_1}, \dots, {\bf u_k}$ が正規直交である」ことを意味しています。この $k$ 本のベクトルを並べたものを $U_k = ({\bf u_1}, \dots, {\bf u_k})$ としておきます。

次に、 ${\bf u_1}, \dots, {\bf u_k}$ に加えて $m - k$ 本のベクトル ${\bf u_{k+1}}, \dots, {\bf u_m}$ をうまくとって正規直交基底を構成します。これはグラム・シュミットの直交化法を用いることで実現することができます。新しく追加した $m - k$ 本のベクトルを並べたものを $U_{m-k} = ({\bf u_{k+1}}, \dots, {\bf u_m})$ とし、 $U = (U_k | U_{m-k}) = ({\bf u_1}, \dots, {\bf u_m})$ とすれば、 $U$ は直交行列になっています。

ここで $U^T A V$ を計算してみましょう。 $U^T A V$ の第 $(i, j)$ 成分を考えると

\begin{eqnarray*}
(U^T A V)_{ij}
&=& (U^T)_{i \cdot} \cdot (A V)_{\cdot j} \\
&=& {\bf u_i}^T {\bf b_j} \\
\end{eqnarray*}

となります。 $j \le k$ のときには

\begin{eqnarray*}
(U^T A V)_{ij}
&=& \sigma_j {\bf u_i}^T {\bf u_j} \\
&=& \sigma_j \delta_{ij}
\end{eqnarray*}

です。 $j > k$ のときには ${\bf b_j} = {\bf 0}$ から $(U^T A V)_{ij} = 0$ が従います。³
よって $U^T A V$ は

\begin{eqnarray*}
U^T A V = 
\left(
\begin{array}{cccc}
\sigma_1 & & & \huge{O} \\
& \ddots & & & \\
\large{O} & & \sigma_k & \\
0 & \dots & \dots & 0 \\
\vdots & & & \vdots \\
0 & \dots & \dots & 0 \\
\end{array}
\right)
\end{eqnarray*}

という $m \times n$ 対角行列となります。これを $\Sigma = U^T A V$ とおけば

$$
A = U \Sigma V^T
$$

という特異値分解を得ることができます。なお、アルゴリズムを見るとわかるように、特異値分解は一意的ではありません。

まとめると、$m \times n$ 行列 $A$ の特異値分解のアルゴリズムは以下のようになります。( $m \ge n$ とします。)

$A^T A$ の固有値分解を実行し、固有値 $\sigma_1^2, \dots, \sigma_n^2$ と対応する固有ベクトル ${\bf v_1}, \dots, {\bf v_n}$ を得る。これらを並べて $V$ とする。
$\sigma_i > 0$ なる $k$ に対して ${\bf u_i} = \frac{1}{\sigma_i} A {\bf v_i}$ とする。(このようにして作成した ${\bf u_i}$ の数を $k$ とおく。)
グラム・シュミットの直交化法によって残り $m - k$ 本の正規直交基底をとり、ステップ2で得たものと合わせて ${\bf u_1}, \dots, {\bf u_m}$ とする。これらを並べて $U$ とする。
対角成分が $\sigma_1^2, \dots, \sigma_n^2$ となる $m \times n$ 対角行列を $\Sigma$ とする。

この手続きを実行することによって、$A$ の特異値分解 $A = U \Sigma V^T$ を得ることができます。ここで $\sigma_i$ の要素を降順に並べる必要はありませんが、各 $i$ について $\sigma_i, {\bf u_i}, {\bf v_i}$ を対応させる必要があることに注意してください。

実装

特異値分解のPythonによる実装はGitHubにあります。

特異値分解

特異値分解を実装するためには固有値分解のメソッドが必要になりますが、ここでは私が実装したヤコビ法(固有値問題) を利用しています。固有値分解のアルゴリズムも非常に重要なので、馴染みのない方はこの機会に勉強しておくと良いと思います。以前ヤコビ法(固有値問題)に関する記事を書いたので、もしよかったら参考にしてください。

後書き

今回は短めの記事になってしまいましたが、このあと特異値分解を利用したアルゴリズムについて記事を書いていく予定です。
特異値分解についてさらに深く理解できるような記事にしていきたいと思います。がんばります。

参考文献・リンク

特異値分解

ダンゴムシでもわかる特異値分解(SVD)
特異値分解を詳しく解説

その他

$m \times n$ 行列 $A$ の対角成分上の要素、すなわち $a_{ii} \ (i = 1, \dots, \min{(m, n)}) $ のことを $A$ の対角成分と呼び、すべての非対角成分が0であるような行列 $A$ を対角行列と呼ぶことにします。 ↩
固有値がすべて非負であるような行列を半正定値行列と言います。詳細は割愛しますが、ある行列 $X$ が存在して $A = X^T X$ と書けるような行列 $A$ は半正定値行列になることが知られています。 ↩
$\| {\bf b_i} \|^2 = \sigma_i^2$ となるため、 $\sigma_i = 0$ の場合にはノルムの公理から ${\bf b_i} = {\bf 0}$ が従います。 ↩