はじめに

前回の記事の続き。
本稿では、みんな大好きウェーブレット変換について記述する。

概要

Waveletsはフーリエ変換の概念をより一般的な直交基底へと拡張し、多重解像度分解を用いることで、部分的に不確定性原理を克服した。さらにウェーブレット変換にも高速計算アルゴリズムが存在し、高次元データへの拡張性を持っている。本稿では、連続ウェーブレット変換を導入し、離散ウェーブレット変換のハイレルな概要をまとめた。

ウェーブレットと多重解像度解析

ウェーブレット解析の基本的な概念は母ウェーブレットと呼ばれる関数$\psi(t)$から始まる。スケール変換の係数を$a$、平行移動の係数を$b$として、次のように母ウェーブレットから異なる一連のウェーブレットを生成することを考える。

\begin{align}
\psi_{a,b}(t)=\frac{1}{\sqrt{a}}\psi\left(\frac{t-b}{a}\right)
\end{align}

この母ウェーブレットの操作は、以下の図(d)のようにセグメントに関数を適用し、それぞれのセグメント内での関数のスケールと平行移動を行うことで、異なる解析手法を実現できる。特に、関数が直交している場合、これらの基底を用いて関数を表現することが可能。

簡単で古典的な例としてHaarウェーブレットを挙げる。

\begin{align}
\psi(t)
=
\begin{cases}
1\quad&\text{for }0\leq t<1/2, \\
-1\quad&\text{for }1/2\leq t<1, \\
0\quad&\text{otherwise.}
\end{cases}
\end{align}

実際にいくつかのパラメータ$a,b$の組み合わせで生成されるウェーブレットを見てみる。

import numpy as np
import matplotlib.pyplot as plt

def haar_wavelet(a, b, t):
    return np.where((b <= t) & (t < b + a / 2), 1,
           np.where((b + a / 2 <= t) & (t < b + a), -1, 0))

# Parameter List
params = [(1, 0), (1/2, 0), (1/2, 1/2)]

# Time axis setting
t = np.linspace(-0.25, 1.25, 1000)

# Create Figure and Axis
fig, axes = plt.subplots(len(params), 1, figsize=(8, 6))
axes[-1].set_xlabel('t')

# Compute and plot Haar wavelets for each parameter set
for i, (a, b) in enumerate(params):
    y = haar_wavelet(a, b, t)
    axes[i].plot(t, y, label=f'Haar Wavelet (a={a}, b={b})')
    axes[i].set_ylabel(f'$\\psi_{{{a},{b}}}$')
    axes[i].grid(True)
    axes[i].legend()

plt.tight_layout()
plt.show()

出力結果は以下の通り。

ここで注目すべきは、各高周波レイヤー(=周波数帯域の信号成分)をその直下のレイヤーの二等分として選択することで、生成されるHaarウェーブレットは直交し、階層的な基底を提供することである。
上の例だと、直下の低周波レイヤー$\psi_{1,0}$の二等分として$\psi_{1/2,0},\ \psi_{1/2,1/2}$を生成したが、$\psi_{1,0}$はこの2つとは直交する。数学的な記述は以下の通り。

\begin{align}
\langle\psi_{a,b},\psi_{a/2,b'}\rangle
=\int_{-\infty}^\infty\psi_{a,b}(t)\psi_{a/2,b'}(t)\,dt
=0
\end{align}

このウェーブレットの直交性は離散ウェーブレット変換（DWT）の拡張へと繋がる。まずは、連続ウェーブレット変換（CWT）を導入する。

\begin{align}
\mathcal{W}_{\psi}(f)(a,b)
=\langle f,\psi_{a,b}\rangle
=\int_{-\infty}^\infty f(t)\overline{\psi}_{a,b}(t)\,dt
\end{align}

これは、関数$\psi(t)$が有界性

\begin{align}
C_{\psi}=\int_{-\infty}^\infty \frac{|\hat{\psi}(\omega)|^2}{|\omega|}\,d\omega
<\infty
\end{align}

を満足する時に有効である。この分数式$\frac{|\hat{\psi}(\omega)|^2}{|\omega|}$は、ウェーブレット関数の周波数成分の寄与度を示し、特に$\frac{1}{|\omega|}$の項は、低周波数領域での寄与が大きくなることを意味する。この積分が有界であることで、ウェーブレット変換が安定となり、信号の情報を適切に捉えることが保証される。

そして、逆変換（iCWT）は以下の通り。

\begin{align}
f(t)=
\cfrac{1}{C_\psi}
\int_{-\infty}^\infty\int_{-\infty}^\infty
\mathcal{W}_{\psi}(f)(a,b)\psi_{a,b}(t)\frac{1}{a^2}\,da\,db
\end{align}

この導出が思いの外、難しかったので以下に残しておく。

まずは、改めてCWTの定義から

\begin{align}
\mathcal{W}_\psi(f)(a,b)
=\int_{-\infty}^\infty f(t)\overline{\psi}_{a,b}(t)\,dt
=f*\overline{\psi}_{a,b}
\end{align}

ここで、$*$は畳み込みを表す。
右辺を$b(t)$などと置いて、着目すると

\begin{align}
b(t)
&=
\cfrac{1}{C_\psi}
\int_{-\infty}^\infty\int_{-\infty}^\infty
\mathcal{W}_{\psi}(f)(a,b)\psi_{a,b}(t)\frac{1}{a^2}\,da\,db \\
&=\cfrac{1}{C_\psi}
\int_{-\infty}^\infty
\mathcal{W}_{\psi}(f)(a,\cdot)*\psi_{a,\cdot}(t)\frac{1}{a^2}\,da \\
&=\cfrac{1}{C_\psi}
\int_{-\infty}^\infty
f*\overline{\psi}_{a,\cdot}*\psi_{a,\cdot}\,\frac{1}{a^2}\,da
\end{align}

$\cdot$は畳み込みが施された変数を示す。これから$b=f$を示していくのだが、上記式を両辺フーリエ変換することでそれを目指す。

\begin{align}
\hat{b}(\omega)
=\cfrac{1}{C_\psi}
\int_{-\infty}^\infty
\hat{f}*\hat{\overline{\psi}}_{a,\cdot}*\hat{\psi}_{a,\cdot}\,\frac{1}{a^2}\,da
\end{align}

非負実数$a$に対する定数倍のフーリエ変換は、$\hat{\psi}_{a,\cdot}(\omega)=\sqrt{a}\hat{\psi}_{1,\cdot}(a\omega)$、そして複素共役については$\hat{\overline{\psi}}(\omega)=\overline{\hat{\psi}(-\omega)}$、さらに$\psi$が実関数であることを用いれば、$\overline{\hat{\psi}(-\omega)}=\overline{\hat{\psi}}(\omega)$なので

\begin{align}
\hat{b}(\omega)
=&\cfrac{1}{C_\psi}
\int_{-\infty}^\infty
\hat{f}*\sqrt{a}\overline{\hat{\psi}}_{1,\cdot}(a\omega)*\sqrt{a}\hat{\psi}_{1,\cdot}(a\omega)\,\frac{1}{a^2}\,da \\
&=\cfrac{\hat{f}}{C_\psi}\int_{-\infty}^\infty
|\overline{\hat{\psi}}_{1,\cdot}(a\omega)|^2\,\frac{1}{a}\,da \\
&=\cfrac{\hat{f}}{C_\psi}\int_{-\infty}^\infty\cfrac{|\hat{\psi}_{1,\cdot}(\xi)|^2}{|\xi|}\,d\xi
=\hat{f}
\end{align}

ウェーブレット$\psi$と有界かつ可積分な関数$\phi$の畳み込み$\psi*\phi$もまたウェーブレットとなる。Haarウェーブレットの他にも有名なウェーブレットは数多く存在し、例えば以下のようなMexican hatウェーブレットなどがある。

\begin{align}
\begin{cases}
\psi(t)&=(1-t^2)e^{-t^2/2},\\
\hat{\psi}(\omega)&=\sqrt{2\pi}\omega^2e^{-\omega^2/2}
\end{cases}
\end{align}

離散ウェーブレット変換

ウェーブレット変換もまた、フーリエ変換/ガボール変換同様に離散変換が存在する。

\begin{align}
\mathcal{W}_\psi(f)(j,k)
=\langle f,\psi_{j,k}\rangle
=\int_{-\infty}^\infty f(t)\overline{\psi}_{j,k}(t)\,dt
\end{align}

ここで、$\psi_{j,k}(t)$は離散ウェーブレット族である。

\begin{align}
\psi_{j,k}(t)=\cfrac{1}{a^j}\psi\left(\frac{t-kb}{a^j}\right)
\end{align}

この離散Haarウェーブレットのように、互いに直交するウェーブレット族を用いれは、それらを基底として関数$f$を表現することが出来る。

\begin{align}
f(t)=
\sum_{j,k=-\infty}^\infty
\langle f(t),\psi_{j,k}(t)\rangle\psi_{j,k}(t)
\end{align}

しかしながら、この明示的な計算は複雑で、それ自体が主題となる論文が多く存在する。ここでは、計算の詳細については言及せず、ハイレベルな"ウェーブレット変換でどんな事が出来るか"、すなわち、信号全体に与えられた形状をスケーリングして変換することにより、時間分解能と周波数分解能の最適なトレードオフを提供する効率的な階層構造で、マルチスケール構造を効率的に抽出することができる、という事をゴールとする。

おわりに

あっけなく終わってしまった。離散ウェーブレットの計算そんなに難しいのか。