ウェーブレット解析についての入門的なメモ

ウェーブレットという言葉ぐらいしか知らなかったので多少調べてメモしておく。

フーリエ解析のおさらい

ウェーブレット解析との比較のためのまとめ。

f(x)とF(ξ)を対称的に解釈することもできるけど、ここではf(x)中心に考える。

フーリエ解析の式は変数や係数の取り方によって細かな違いが出てくるけど、以下の表はその一例。

xのとる範囲 関数f(x)に、周波数ξの成分がどれだけ含まれているか(=F(ξ))を求める 関数f(x)を、周波数成分の重ね合わせで表現する
離散フーリエ変換  x \in \{0, \cdots , N-1 \}  F(\xi) = \frac{1}{N} \sum_{x = 0}^{N-1} f(x) e^{ - 2\pi i \frac{\xi x}{N}}  f(x) = \sum_{\xi = 0}^{N-1} F(x) e^{  2\pi i \frac{\xi x}{N}}
離散時間フーリエ変換  x \in Z  F(\xi) = \sum_{x=- \infty}^{\infty} f(x) e^{-2 \pi i \xi x}   f(x) = \int\limits_0^1 F(\xi) e^{2 \pi i \xi x} d\xi
フーリエ展開  0 \le x \lt 1  F(\xi) = \int\limits_0^1 f(x) e^{-2 \pi i \xi x} dx  f(x) = \sum_{\xi=- \infty}^{\infty} F(\xi) e^{2 \pi i \xi x}
フーリエ変換  x \in R  F(\xi) = \int\limits_{-\infty}^{\infty} f(x) e^{-2 \pi i \xi x} dx  f(x) = \int\limits_{-\infty}^{\infty} F(\xi) e^{2 \pi i \xi x} d\xi

ウェーブレット解析と比較しやすいように \psi(x) = e^{2 \pi i x} と置くと次のようになる。

離散フーリエ変換  F(\xi) = \frac{1}{N} \sum_{x = 0}^{N-1} f(x) \bar{\psi(\frac{\xi x}{N}})  f(x) = \sum_{\xi = 0}^{N-1} F(x) \psi  (\frac{\xi x}{N})
離散時間フーリエ変換  F(\xi) = \sum_{x=- \infty}^{\infty} f(x) \bar{\psi(\xi x)}   f(x) = \int\limits_0^1 F(\xi) \psi(\xi x) d\xi
フーリエ展開  F(\xi) = \int\limits_0^1 f(x) \bar{\psi( \xi x)} dx  f(x) = \sum_{\xi=- \infty}^{\infty} F(\xi) \psi(\xi x)
フーリエ変換  F(\xi) = \int\limits_{-\infty}^{\infty} f(x) \bar{\psi(\xi x)} dx  f(x) = \int\limits_{-\infty}^{\infty} F(\xi) \psi(\xi x) d\xi

これを見ると、変換のための関数ψ(x)(指数関数、三角関数)に対して、その引数を定数倍する(=関数の形を横方向に縮小する)という作業をしていることが判る。ψ(ξx)は形が圧縮される(ξが大きくなる)につれて、より高周波の波になっていく。この波をf(x)をかけあわせて積分(or和)を取ることで、f(x)の特定の周波数成分だけをより分けている。

ウェーブレット変換の式

一方、ウェーブレット変換の式は次のようになっている。

離散ウェーブレット変換  F(k,l) =\int_{-\infty}^{\infty} f(x) \bar{ \frac{1}{\sqrt{A^k}} \psi(\frac{x}{A^k}- B l)} dx  f(x) = \sum_{k,l} F(k,l) \frac{1}{\sqrt{A^k}} \psi(\frac{x}{A^k}- B l)
(ψ/√が正規直交基底になっている場合)
連続ウェーブレット変換  F(a,b) = \int_{-\infty}^{\infty} f(x) \bar{\frac{1}{\sqrt{|a|}} \psi(\frac{x-b}{a})}dx  f(x) = \frac{1}{C_{\psi}}\int \int F(a,b) {\frac{1}{\sqrt{|a|}}\psi(\frac{x-b}{a}) db \frac{da}{a^2} }
ここでのψ(x)は、指数関数・三角関数ではなくて別の何らかの関数。一番大きな違いは、フーリエ変換のψ(x)は絶対値|ψ(x)|が常に一定(三角関数の場合はcosとsinを組にして考える)なのに対して、ウェーブレット変換のψ(x)は遠くの方では0に近づく関数だということ。また指数関数・三角関数とは違って、ψは直交関係を作るような関数である必要はなくて、もっとゆるい条件を満たすだけでよい。
ただしそのゆるい条件から導かれる \textstyle{\int_{-\infty}^{\infty} \psi(x) dx = 0}という条件は重要。この条件から

  • ψ(x)は、局在化した(=遠くでは0に近づく)関数である(そうでないと条件の積分の値が確定しないから)。
  • ψ(x)は、おおざっぱに見れば波のような形になる(正の部分の面積と負の部分の面積が等しくなるため。ただし波形は連続とは限らない)。

となることがわかる。つまりψ(x)は、局在化した波のような形の関数でないといけない(ウェーブレット関数と呼ばれる)。
連続ウェーブレット変換の場合、条件を満たすψを使ってF(a,b)を求めれば逆変換をすることができる。
一方、離散ウェーブレット変換の場合、うまいψを取って、さらにうまく離散点を取らないと元の関数f(x)が再現できないかもしれない。あるいは再現はできるけど必要以上にたくさんの点のデータを取った冗長な形なっているかもしれない。おそらく一番都合が良いのは、 {\frac{1}{\sqrt{A^k}} \psi(\frac{x}{A^k}- B l)\quad (k,l \in Z)}が正規直交基底になるようにψをとることで、この場合はフーリエ展開と同じような逆変換の式になる。
連続、離散いずれにしてもなるべく都合の良いψを取ることが重要で、ψの選び方や作り方の話がウェーブレット解析の大きな部分になっている。おそらくは。

ウェーブレット変換の値F(a,b)は何を表しているか

ウェーブレット関数ψ(x)が x=b_0の周りに局在化した関数だとすると、 \psi(\frac{x-b}{a})は、 x= b_0+bに局在化している。よってf(x)の値のうち x=b_0+bの周りのものだけが {F(a,b) = \int f(x) \bar{\frac{1}{\sqrt{|a|}} \psi(\frac{x-b}{a})}dx }に影響を与える。なのでF(a,b)は、 x= b_0+bの周りでのf(x)の何らかの特徴を見ていることになる。
一方aはフーリエ変換から類推すると、周波数に関する量を表していると予想できる。ただしψ(x)は三角関数ではなくさまざまな周波数の影響を合わせ持った関数 \textstyle{\psi(x) = \int \hat{\psi}(\xi) e^{2 \pi i \xi x}d\xi}なので、正確な周波数成分には対応しない。
それでもψ(x)のフーリエ変換 \hat{\psi}(\xi)(これはψ(x)の周波数成分をあらわしている)のピークが \xi_0にあるとすると \hat{\psi(\frac{x-b}{a})}(\xi)のピークは \frac{\xi_0}{a}にあるので、F(a,b)の値はf(x)のうち周波数 \frac{\xi_0}{a}の周辺の特徴をとらえていると考えることができる。
よってウェーブレット変換で得られるF(a,b)は、関数f(x)のうち

  • 位置x=b_0+bの周り
  • 周波数 \xi = \frac{\xi_0}{a}の周り

の特徴から決まる(あらかじめψを平行移動・拡大縮小しておけば「位置b、周波数1/aの周り」になる)。ただしaが小さくなる(周波数が大きくなる)ほど、見ている空間領域は狭くなり見ている周波数領域は広くなる(aが大きくなると \psi(\frac{x-b}{a})は空間領域での幅が狭くなり、周波数領域での幅が広がるので。高周波成分≒急激な変化なので空間の狭い範囲だけを見れば良い。低周波成分≒ゆるやかな変化なのでより広い範囲を見る必要がある)。
もちろんψ(x)の取り方によってどれくらいの範囲をどのように見るかが変わりF(a,b)の値も変わるけれど、定性的にはF(a,b)の値(の絶対値)が大きければ大きいほどf(x)は位置bで周波数1/aの成分をたくさん含んでいると言える。たぶん。

多重解像度解析

ウェーブレット解析についての文章を流し見ると「画像が4分割されていて、4つの分割のひとつは元画像を縮小したもので残りの部分はよく分からないデータになっているもの」を見かけることがある(縮小された画像はさらに再分割・縮小化されていたりもする)。たとえばWikipediaWavelet transformに載っている画像。あれは何をしているかについて。
直交基底を使った離散ウェーブレット変換は、元の関数f(x)を
 f(x) = \sum_{k} \sum_{l}F(k,l) \left(\frac{1}{\sqrt{A^k}} \psi(\frac{x}{A^k}- B l) \right)
のように、周波数・空間成分に分解する。このときaとbの値は (a,b)=(A^k, A^{-k}B l)で、高周波(aが小さい)になるほど空間成分の分解は細かくなる(データの数が増える)。またフーリエ展開との類推によると、低周波成分の方が相対的に重要なデータになる。
それを考えるとF(k,l)の値は、低周波成分から順番に計算したくなる。でも基底関数がうまい性質を持っている場合、逆に高周波成分から計算した方が効率良く計算できる。
まずf(x)に含まれている一番高周波の成分をk=1としておく。そのとき、もっとも高周波の成分\Psi_1(x)とそれ以外の残り \Phi_1(x)に分解する。
 f(x)  = \Phi_1(x) + \Psi_1(x) = \sum_{k>1} \sum_{l}F(k,l) \left(\frac{1}{\sqrt{A^k}} \psi(\frac{x}{A^k}- B l) \right) + F(1,l) \left(\frac{1}{\sqrt{A^1}} \psi(\frac{x}{A^1}- B l) \right)
ここで高周波成分を取り除くという作業は、データを少し平均化して細かな変化を取り除いていると解釈できる。つまり f(x)=\Phi_1(x)+\Psi_1(x)は、「平均化したデータ+細かい変化」に分解している。
さらにウェーブレット変換に使っている基底関数ψ(x)がうまい性質を持っている場合は \Phi_1(x)も(基底関数とは別の関数φ(x)を使って)成分表示することができるため効率よく計算ができる。そして同様の分解をさらに繰り返せば変換が終了する。ウェーブレット変換の話で見かける4分割化された図は、2次元画像でこの分解をしている(x軸、y軸のそれぞれで分解するから4分割になる)。
このような分解のやり方を多重解像度解析と呼ぶ。これが可能な場合、背後に特別な空間構造が存在しているので、この考え方は性質の良い基底を構成するためにも利用される。

逆変換の式

連続ウェーブレット変換の逆変換の式 f(x) = \frac{1}{C_{\psi}}\int \int F(a,b) {\frac{1}{\sqrt{|a|}}\psi(\frac{x-b}{a}) db \frac{da}{a^2} }はいくらかややこしい式になっていて、これだけだと離散変換の式 f(x) = \sum_{k,l} F(k,l) \frac{1}{\sqrt{A^k}} \psi(\frac{x}{A^k}- B l)  との対応がよく判らない。

離散版の逆変換式から連続版の逆変換式を導く

離散版ではF(a,b)のうち飛び飛びの点だけを使っているので、この点をどんどん細かくしていけば連続の式に移行するだろうという方針が最初に思いつく。でも使っている点は (a,b)=(A^k, A^{-k}B l)でAとBは固定なので、分割を細かくできない(AやBを変えると正規直交基底を成すという性質が崩れてしまう)。
どうするかというと、 \frac{1}{\sqrt{A^k}} \psi(\frac{x}{A^k}- B l)が正規直交基底になっているならば、任意の定数c、dについて \frac{1}{\sqrt{cA^k}} \psi(\frac{x-d}{cA^k}- B l )も正規直交基底になることに注目する。この新しい基底で変換をおこなった場合、 (a,b)=(cA^k, cA^{-k}B l+d)という元の点とはズレた点のF(a,b)を計算することになる。この点を元の点と合わせれば、分割の間隔を小さくすることができる。
自然数Nを取り c=a^{\frac{m}{N}},d=\frac{n}{N}\, (m,n=0,\ldots,N-1)のそれぞれについてF(a,b)を求めれば、f(x)を表す式がN2個できる。それらを足してN2で割れば
 f(x) = \frac{1}{N^2}\sum F(a,b) \frac{1}{\sqrt{a}}\psi(\frac{x-b}{a})
が得られる。この式での(a,b)の分割の間隔は (\Delta a, \Delta b) = ((A^{\frac{1}{N}}-1)a, \frac{aB}{N})になるので、
 f(x) = \sum F(a,b) \frac{1}{\sqrt{a}}\psi(\frac{x-b}{a})  \frac{\Delta a \Delta b}{a^2 B N(A^{\frac{1}{N}}-1)}
と書けて、N→∞で間隔を狭めていけば
 f(x) = \frac{1}{B \log A} \int \int F(a,b) \frac{1}{\sqrt{a}}\psi(\frac{x-b}{a})  \frac{d a }{a^2}db \quad =\frac{1}{C} \int \int F(a,b) \frac{1}{\sqrt{a}}\psi(\frac{x-b}{a})  \frac{d a }{a^2}db
と連続の場合の逆変換式になる(ここではaの積分範囲は正になっているけど、実数全体で考えることもできる)。

連続版の逆変換式

上でおこなった離散から連続へ変形は、離散版の式が成り立つこと、つまりψから正規直交基底が作れることを前提にしている。でも連続版の式自体はそのような前提なしでも成り立つ。
f(x)、ψ(x)のフーリエ変換 \hat{f}(\xi) \hat{\psi}(\xi)を使って変形していくと
 F(a,b) = \int \hat{f}(\eta) \sqrt{|a|}\bar{\hat{\psi}(a\eta)} e^{2\pi i b \eta} d \eta

 F(a,b) \frac{1}{\sqrt{|a|}}\psi(\frac{x-b}{a}) = \int \int \hat{f}(\eta)|a| \hat{\psi}(a \mu) \bar{\hat{\psi}(a \eta)} e^{2 \pi i \mu x} e^{2 \pi b(\eta - \mu)}d\eta d\mu

 \int F(a,b) \frac{1}{\sqrt{|a|}}\psi(\frac{x-b}{a})db =  \int \hat{f}(\eta)|a| |\hat{\psi}(a \eta)|^2 e^{2 \pi i \eta x} d\eta

 \int \int F(a,b) \frac{1}{\sqrt{|a|}}\psi(\frac{x-b}{a})db \frac{da}{a^2}=  \int \hat{f}(\lambda )  e^{2 \pi i \lambda x} d\lambda   \cdot \int  \frac{|\hat{\psi}( \xi)|^2}{|\xi|} d\xi
となり、 \textstyle{\int\frac{|\hat{\psi}( \xi)|^2}{|\xi|} d\xi}=C<\inftyの場合は逆変換の式が成り立つ。