LoginSignup
4
3

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

More than 1 year has passed since last update.

@arcadia13(arcadia 13)

MATLABで特異値分解を理解する。

Last updated at Posted at 2022-08-10

10回目のMATLAB記事です。

今回のメイン「特異値分解（Singular Value Decomposition）」です。

特異値分解は、固有値分解の親戚みたいな存在です。

固有値分解は正方行列に対して定義されました。

じゃあ非正方行列はどうするのか、と疑問に思いますよね。

それが、今回説明する特異値分解です1。

本記事では、MATLABでの特異値分解の実装と、その応用例である画像の低ランク近似を実装します。

本記事はある程度の線形代数の知識を前提にしています。

特異値分解

特異値と特異値分解

ある行列 $A \in \mathbb{C}^{m\times n}(m \neq n)$ を考えます。

このとき

  • 対角成分が非負の実数、非対角成分が $0$ の行列 $\Sigma \in \mathbb{R}^{m\times n}$
  • 二つのユニタリ行列 $U \in \mathbb{C}^{m\times m},\ V \in \mathbb{C}^{n\times n}$

を用いて行列 $A$ を

\begin{equation}
  A=U\Sigma V^{*}
\end{equation}

と分解することができます。

ただし、$V^{*}$ は $V$ の共役転置を表します2。

これを行列 $A$ の特異値分解といい、$\Sigma$ の対角成分の中で $0$ でないものを特異値といいます。

また、もし $A$ が実行列の場合、$U,V$ が実直交行列となり、特異値分解は

\begin{equation}
  A=U\Sigma V^{T}
\end{equation}

となります。

共役転置が単純に転置になっただけです。

特異値分解の具体例

先に具体例を示します。

以下の行列の特異値分解を考えます。

\begin{equation}
A=\begin{pmatrix}
2 &  2 & 2 & 2\\
1 & -1 & 1 & -1\\
-1 & 1 & -1 & 1
\end{pmatrix}
\end{equation}

この行列の特異値分解は以下の通りです。

\begin{equation}
A=
\begin{pmatrix}
1 &  0 & 0\\
0 & \frac{1}{\sqrt{2}} & \frac{1}{\sqrt{2}}\\
0 & \frac{-1}{\sqrt{2}} & \frac{1}{\sqrt{2}}
\end{pmatrix}
\begin{pmatrix}
4 & 0 & 0 & 0\\
0 & 2\sqrt{2} & 0& 0\\
0 & 0 & 0 & 0
\end{pmatrix}
\begin{pmatrix}
\frac{1}{2} & \frac{1}{2} & \frac{1}{2} & \frac{1}{2}\\
\frac{1}{2} & \frac{-1}{2} & \frac{1}{2} & \frac{-1}{2}\\
\frac{1}{2} & \frac{1}{2} & \frac{-1}{2} & \frac{-1}{2}\\
\frac{1}{2} & \frac{-1}{2} & \frac{-1}{2} & \frac{1}{2}
\end{pmatrix}
\end{equation}

特異値は $4, 2\sqrt{2}$ です。

特異値の性質

特異値には、$A^{*}A$ または $AA^{*}$ の固有値の平方根に等しいという性質があります。

つまり、特異値を手計算する時は、$A^{*}A$ または $AA^{*}$ の固有値を求めて、それぞれの平方根を計算すればいいのです。

特異値分解の別表現

ここからは、$A$ が実行列だとして話を進めます3。

実直交行列 $U,V$ の各列をそれぞれ

\begin{equation}
\boldsymbol{u}_j \in \mathbb{R}^{m}, \boldsymbol{v}_k \in \mathbb{R}^{n}\ (j=1,2,\dots,m,\ k=1,2,\dots,n)
\end{equation}

と列ベクトルで表します。

すると、$U,V$ は

\begin{equation}
U=\begin{pmatrix}
\boldsymbol{u}_1,\ \boldsymbol{u}_2,\ \dots\ ,\ \boldsymbol{u}_m 
\end{pmatrix}\\

V=\begin{pmatrix}
\boldsymbol{v}_1,\ \boldsymbol{v}_2,\ \dots\ ,\ \boldsymbol{v}_n 
\end{pmatrix}
\end{equation}

となります。

特異値分解で $V^T$ が出てきますが、これは

\begin{equation}
V^{T}=
\begin{pmatrix}
\boldsymbol{v}_{1}^T \\
\boldsymbol{v}_{2}^T \\
\vdots\\
\boldsymbol{v}_{n}^T
\end{pmatrix}
\end{equation}

のように、列ベクトルの転置（すなわち行ベクトル）を下にどんどん並べていったものです。

$\boldsymbol{u}_j, \boldsymbol{v}_k$ はそれぞれ、左特異ベクトルと右特異ベクトルといいます。

名前の単純な覚え方としては、$\Sigma$ の左側、右側にきているからと覚えてください。

さて、この二種類の特異ベクトルを使って、行列 $A$ の特異値分解は、以下のように書き換えることができます。

ただし、$\sigma_i$ は $\Sigma$ の対角成分です。

\begin{equation}
A=\sum_{i=1}^{R}\sigma_{i}\boldsymbol{u}_{i}\boldsymbol{v}_{i}^{T},\ R=\min(m,n)
\end{equation}

ピンとこない方は、先ほどの具体例に対してこの書き換えを行ってみてください。

確かにこのような式で表せることが分かると思います。

また、右辺の項

\sigma_{i}\boldsymbol{u}_{i}\boldsymbol{v}_{i}^{T}

は、行列 $A$ と同じサイズの「行列」です。

一見、そうは見えないかもしれませんが、左から

  • $\sigma_i$：スカラ
  • $\boldsymbol{u}_{i}$：列ベクトル
  • $\boldsymbol{v}_{i}^{T}$：列ベクトルの転置（すなわち行ベクトル）

の積なので、行列です。

特異値分解は行列の和として表すことができるのです。

行列の低ランク近似

行列 $A$ を、それよりも低ランクの行列で近似することを考えます。

具体的には、特異値分解の別表現

\begin{equation}
A=\sum_{i=1}^{R}\sigma_{i}\boldsymbol{u}_{i}\boldsymbol{v}_{i}^{T},\ R=\min(m,n)
\end{equation}

において、$r<R$ を満たす $r$ で打ち切ってしまう、というものです。

つまり、行列の足し算を適当なところで打ち切るのです。

そして出来上がった行列は、完全な元の行列ではなく、それに近しい行列にすぎません。

ゆえに近似となります。

これを、行列 $A$ の低ランク近似といいます。

特異値分解の実装

固有値分解と同様、特異値分解の実装はすごく簡単です。

svd() コマンドを使うだけです。

sample1.m
A = [2,2,2,2; 1,-1,1,-1; -1,1,-1,1];

% 特異値分解
[U,S,V] = svd(A);

U
S
V
U*S*V' %Aの復元
実行結果
U =

   -1.0000         0         0
         0    0.7071    0.7071
         0   -0.7071    0.7071


S =

    4.0000         0         0         0
         0    2.8284         0         0
         0         0    0.0000         0


V =

   -0.5000    0.5000    0.1861    0.6822
   -0.5000   -0.5000    0.6822   -0.1861
   -0.5000    0.5000   -0.1861   -0.6822
   -0.5000   -0.5000   -0.6822    0.1861


ans =

    2.0000    2.0000    2.0000    2.0000
    1.0000   -1.0000    1.0000   -1.0000
   -1.0000    1.0000   -1.0000    1.0000

分かりやすいことに、戻り値を受け取る変数 U、S、V は左から $U,\Sigma,V$ となっています。

また、これらから元の行列 $A$ を復元できます。

低ランク近似の実装

低ランク近似は、以下のように実装します。

sample2.m
A = [2,2,2,2; 1,-1,1,-1; -1,1,-1,1];

% 特異値分解
[U,S,V] = svd(A);

% 低ランク近似
r = 1; % ここを1～3に自由に変更してください。

A_cal = 0;
for i=1:r
    A_cal = A_cal + U(:,i)*S(i,i)*V(:,i)';
end
A_cal
実行結果
A_cal =

     2     2     2     2
     0     0     0     0
     0     0     0     0

単純なfor文によって低ランク近似を実装できます。

for文の中身が何をしているのかというと、変数 A_cal に

\sigma_{i}\boldsymbol{u}_{i}\boldsymbol{v}_{i}^{T}

を次々に足し上げているのです。

そして最終的に、A_cal は行列 $A$ の低ランク近似になります。

また、変数 r を1、2に変えることで、低ランク近似を行うことができます。

r を3にすれば元の特異値分解です。

試しに r を3にして実行してみました。

実行結果
A_cal =

    2.0000    2.0000    2.0000    2.0000
    1.0000   -1.0000    1.0000   -1.0000
   -1.0000    1.0000   -1.0000    1.0000

確かに元に戻っています。

応用例：画像の低ランク近似

特異値分解も、固有値分解と同様に様々な応用例があります。

今回はその中で、低ランク近似の応用例である「画像の低ランク近似」を紹介します。

画像は行列である

画像をコンピュータで読み込む時、それは行列として読み込みます。

例えばグレー画像の場合、各要素が0から255の値になり、その値は各ピクセルの白黒度合いを表します。

我々人間にとって、画像とは絵や写真という認識でしょう。

一方で、コンピュータにとって画像は行列、つまり数字を並べたものという認識なのです。

また、画像は非正方行列であることが多く、そのため特異値分解を前提にした方が都合がいいのです。

題材にした画像

今回は以下の、猫のフリー画像を題材にします。

これはRGB画像なので、特異値分解を行いやすいようにグレー画像へ変更します。

全くの蛇足ですが、私はドが付くほどの愛猫家です（笑）。

低ランク近似の実装

低ランク近似は以下のように実装します。

言うまでもないかもしれませんが、MATLABで簡単な画像処理を実装することになります。

sample3.m
clear; close all;

% RGB画像の読み込み
% この sample3.m と画像データを同じフォルダに入れてください
im_rgb = imread('cat.jpg');

% 元画像の表示
figure(1)
imshow(im_rgb)

% RGB画像をグレー画像へ変換
im_gr = rgb2gray(im_rgb);

% グレー画像の表示
figure(2)
imshow(im_gr)

% unit8配列をdouble配列に変更
% unit8形式だと特異値分解する上で都合が悪いのでdouble形式に変換
A = im2double(im_gr);

% 特異値分解
[U,S,V] = svd(A);

% 低ランク近似
r = 90; % ここを1から333まで自由に変更してください。
A_cal = 0;
for i=1:r
    A_cal = A_cal + U(:,i)*S(i,i)*V(:,i)';
end

% 低ランク近似した画像
figure(3)
imshow(A_cal)

% 低ランク近似第1項
figure(4)
imshow(U(:,1)*S(1,1)*V(:,1)')

% 低ランク近似第2項
figure(5)
imshow(U(:,2)*S(2,2)*V(:,2)')

実行すると、五枚の画像が表示されると思います（Figure1からFigure5）。

処理の順番の都合上、一枚目が一番下、二枚絵がその一つ上......となっています。

一枚目の画像は元の画像（RGB画像）表示しただけですので割愛します。

二枚目と三枚目の画像

二枚目は元の画像をグレー画像にしたもので、三枚目が低ランク近似によりそれを近似した画像です。

上が二枚目、下が三枚目です。

若干、三枚目の画質が荒いもののパッと見、違いは分かりません。

もちろん、r の値を上げていけばどんどん再現度が上がっていきます。

実用上、パッと見、違いが分からないのであればそれを無視することがよくあります。

なお、r の値は333まで上げられます。

333の根拠は、画像の行数（$m$）と列数（$n$）のうち、小さい方がこの値だったからです。

四枚目と五枚目の画像

次に四枚目と五枚目。

「何だこれは？」と思うでしょう。

特異値分解は行列の和として表せると説明しました。

また、画像はコンピュータにとって行列です。

つまり、元のグレー画像は複数のグレー画像の和によって表すことができるのです。

四枚目と五枚目は、その第一項と第二項の画像を表したものです。

これらをどんどん足し上げることで、元のグレー画像を復元できるのです。

まとめ

今回は、特異値分解の実装と、その応用例である画像の低ランク近似をMATLABで実装しました。

私自身、画像処理に詳しくないため、大規模なもの、複雑なものは実装しませんでした。

ただ、数式をそのまま読むよりも特異値分解を直感的に理解できるのではないかと思います。

また、特異値分解は主成分分析という手法で使われます。

主成分分析は統計学や機械学習で重要な手法で、興味のある方は調べてみてください。

分解ではなく、特異値であれば、制御理論のロバスト制御という制御方式があって、そこで活用しています。

めっちゃ難しいですが、これも興味があれば調べてみてください。

終わりに

「特異」は、日本では良い意味で使われることの少ない言葉です。

この言葉の意味を調べてみると

  1. 特別に他とちがっていること。また、そのさま。
  2. 特にすぐれていること。また、そのさま。

コトバンク：デジタル大辞泉「特異」の解説より

だそうです。

後者の意味を採用すれば、とても良い言葉だと感じます。

前者の意味であっても、別に悪いニュアンスを持ってるわけではありません。

なぜなら、単に「異なっている」という状態を表しているにすぎないからです。

平均して、プラスのイメージがある言葉ではないでしょうか。

では、「異」という漢字を見た時にどんな言葉を連想しますか？

おそらく「異なる」や「異常」が、大体を占めるのではないでしょうか。

日本では「異常」は、良い言葉ではありません。

元々、状態異常のように何かしら悪い状態を表す言葉ですからね。

「特異」は「異常」のイメージに引っ張られて、良い意味で使われることが少ないのかな、と私は思います。

次に疑問に思うのは特異値という名前の由来です。

「特異な値ってなんだよ？」と疑問になりませんか？

おそらく、「特」から「特別な値」という意味にしたいのでしょう。

ただ、何が特別なのか、辞書的に考えれば何が優れているのか、異なるのか。

私はいまいちピンときていません。

まだまだ私は勉強不足ですね。

ぜひ名付け親に真意を訊きたいところです。

よろしければ次回の記事も読んでくださると大変嬉しいです。

※本記事に対する改善点や修正点、またはこんな事が知りたいといったご意見がありましたらぜひご連絡ください。

参考HP

  1. 学研ホールディングス, 高校数学の美しい物語, 特異値分解の定義，性質，具体例, 2021/03/07更新, 2022/8/9閲覧
    https://manabitimes.jp/math/1280

  2. 機械学習ナビ, 【線形代数】特異値分解とは?例題付きで分かりやすく解説!!, 森田 大樹, 2022/03/14更新, 2022/8/5閲覧
    https://nisshingeppo.com/ai/singular-value-decomposition/

  3. Qiita, 特異値分解を詳しく解説, @kidaufo, 2019/02/26投稿, 2021/01/03更新, 2022/8/6閲覧
    https://qiita.com/kidaufo/items/0f3da4ca4e19dc0e987e

  4. Think IT, 分解すると見える世界 ー特異値分解ー, 森田 大樹, 2019/10/23投稿, 2022/8/5閲覧
    https://thinkit.co.jp/article/16884

  5. Qiita, 特異値分解による画像の低ランク近似, @kaityo256, 2016/11/17投稿, 2020/06/03更新, 2022/8/7閲覧
    https://qiita.com/kaityo256/items/48de63526b469235d16a

  1. 実は、特異値分解は正方行列に対しても行うことができます。 ↩

  2. 共役転置は複素行列に対して行う演算で、各要素の共役複素数をとった後、さらに転置するというものです。実行列の転置に相当する演算で、実際、実行列を共役転置した行列は元の行列の転置行列と一致します。 ↩

  3. 複素行列だとしても、実直交行列をユニタリ行列、転置を共役転置に置き換えるだけで同様の議論が成立します。 ↩

4
3
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
Sign upLogin
4
3

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?