はじめに
テンソルを扱っていると、「あれ？いまどの足がどこにあるんだっけ？」とごちゃごちゃになることが多い。というわけで、HOSVDによる画像の低ランク近似を題材に、一つ一つ確認しながらステップを進めてみる。
以下、実習を想定し、過程をしつこいくらい詳細に書く。
ステップ・バイ・ステップ
ステップ1 画像の読み込みまで
以下、表示の都合上Pythonの対話モードで書いているが、実際にはipythonでやったほうがいろいろ便利。
まずは必要となるライブラリをインポートしよう。
>>> import numpy as np
>>> from scipy import linalg
>>> from PIL import Image
次に適当な画像を読み込む。ちょうど手元に手頃なサイズの画像があったので、それを使うことにしよう。
これはItanum2の精霊馬である。
これを読み込んで、サイズを確認してみよう。
>>> img = Image.open("uma.jpg")
>>> img.size
(480, 340)
480X340ピクセルの画像データだ。後のため、高さと幅を変数に落としておく。
>>> w, h = img.size
さて、イメージデータをNumpyの配列にしてみる。
>>> X = np.asarray(img)
>>> X.shape
(340, 480, 3)
順番が「高さ、幅、3(RGB)」となっていることに注意。
ステップ2 足を入れ替えてSVDの準備
これからSVDをするのだが、SVDをする足を一番右に持ってきて、さらに3階のテンソルを行列の形にする必要がある。足の入れ替えにはtransposeを使う。
再度Xの次元を確認する。
>>> X.shape
(340, 480, 3)
これはいま「高さ、幅、色」となっているが、まず幅について足を潰したいので、「高さ、色、幅」の順番にしなければならない。「高さ、幅、色」の順番が(0,1,2)なので、目標とする順番は(0,2,1)である。
>>> X.transpose(0,2,1).shape
(340, 3, 480)
順番が入れ替わった。これをさらに、左２つの足をまとめて行列にする。
>>> X.transpose(0,2,1).reshape(h*3,w).shape
(1020, 480)
無事に、(高さと色, 幅)の次元を持つ行列になった。これをX1として保存しておこう。
>>> X1 = X.transpose(0,2,1).reshape(h*3,w)
同様に、一番右に高さの次元を持つ行列を作りたい。もう一度Xの次元を見る。
>>> X.shape
(340, 480, 3)
これを(480, 3, 340)にしたい1。そのためには(1,2,0)の順番にすれば良い。
>>> X.transpose(1,2,0).shape
(480, 3, 340)
同様に、左の足を2本まとめて行列にして、これをX2に保存しておこう。
>>> X2 = X.transpose(1,2,0).reshape(w*3,h)
>>> X2.shape
(1440, 340)
右の足の次元が高さになっていますね。ここまででSVDの準備完了。
ステップ3 SVDしてコアテンソルを作る
まずX1をSVDする。さらにそれぞれのサイズを確認してみる。
>>> U, s, A1 = linalg.svd(X1)
>>> U.shape
(1020, 1020)
>>> s.shape
(480,)
>>> A1.shape
(480, 480)
さて、A1の情報を5%に圧縮する。もともと480行だったので、上位24行だけ取ってくる。これをa1としよう。
>>> a1 = A1[:24, :]
>>> a1.shape
(24, 480)
同様にX2をSVDして、A2を上位17行(340の5%)だけとってくる。
>>> U, s, A2 = linalg.svd(X2)
>>> A2.shape
(340, 340)
>>> a2 = A2[:17,:]
>>> a2.shape
(17, 340)
ここでえられたa1, a2の転置を取ると、幅、高さをつぶす圧縮行列になる。
>>> a1.T.shape
(480, 24)
>>> a2.T.shape
(340, 17)
それぞれ480→24、340→17に潰す行列となっていることがわかる。これをXの適当な足と内積を取ればコアテンソルを得る。
まず、Xにa1.Tをかける。どの足とどの足をつぶせば良いのかは、shapeを見ればわかる。
>>> X.shape
(340, 480, 3)
>>> a1.T.shape
(480, 24)
一番左を0番目として、Xの1番目とa1.Tの0番目の足を潰せば良い。これをX3に代入する。
>>> X3 = np.tensordot(X, a1.T, (1,0))
>>> X3.shape
(340, 3, 24)
無事に480→24に圧縮された。ここで注意すべきは、tensordot後の足の順番。慣れるまでは新しいテンソルのどの足が何になっているかを毎回shapeで確認したほうが良いと思う。
次に、X3の高さを表す足にa2.Tをかけて圧縮する。あらためてどれとどれをつぶすべきか調べる。
>>> X3.shape
(340, 3, 24)
>>> a2.T.shape
(340, 17)
足の大きさが全部違う場合は、同じ次元のところをつぶせばよいとすぐわかって便利。この場合は0番同士をつぶせば良いので、
>>> x = np.tensordot(X3,a2.T,(0,0))
>>> x.shape
(3, 24, 17)
無事に、Xの次元を削減したコアテンソルxが得られた。
最終的に、Xから、幅の圧縮/復元行列a1, 高さの圧縮/復元行列a2, コアテンソルxに圧縮された。この圧縮率を調べてみよう。
>>> float(x.size + a1.size + a2.size)/X.size
0.03783496732026144
ということで、もともとのデータの3.8%くらいまで圧縮できたことになる。
ステップ4 情報の復元
いま、手元にはコアテンソルxと、幅、高さの復元行列a1,a2がある。コアテンソルにa1,a2をかければ、元のテンソルの近似テンソルYが得られるので、それを作って見よう。
まずは次元の確認。
>>> x.shape
(3, 24, 17)
>>> a1.shape
(24, 480)
xの1番とa1の0番をつぶせばよいから、
>>> x2 = np.tensordot(x,a1,(1,0))
>>> x2.shape
(3, 17, 480)
幅が24から480に復元された。同様に高さを復元したものをx3とする。
>>> x3 = np.tensordot(x2,a2,(1,0))
>>> x3.shape
(3, 480, 340)
ただし、このままでは足の順番が違う。もう一度Xを見てみる。
>>> X.shape
(340, 480, 3)
Xと見比べれば、足を(2,1,0)の順番にすればよいことがわかるので、
>>> Y = x3.transpose(2,1,0)
>>> Y.shape
(340, 480, 3)
これが、コアテンソルxから、復元行列a1, a2を用いて復元されたテンソルである。
これで元のデータと同じ次元になったので、Imageオブジェクトを作ることができる。ただし、Image.fromarrayにつっこむ前にuint8型にしておく必要がある。
>>> img2 = Image.fromarray(np.uint8(Y))
>>> img2.save("approximated.jpg")
こうしてできた復元画像がこんな感じ。
まぁ、「元が画像である」という情報をまったく使わずに3.8%まで圧縮したデータから復元したわりにはそこそこ復元できたかな、という感じ。
まとめ
HOSVD(Higher Order Singular Value Decomposition)を使って画像を圧縮、復元してみた。IPython使うと、タブ補完とか効いてすごく便利。テンソルの足を潰したあと、何番目が何を表す足だったか混乱しがちだが、いちいちshapeで次元を調べながら作業すると間違い/混乱も減る気がする。
しっかしPythonは便利だね・・・2。