多重解像度解析

多重解像度解析（たじゅうかいぞうどかいせき、）とは、2倍毎の解像度のウェーブレットを用いて離散ウェーブレット変換により解析する手法。スケーリング関数で基底展開された信号列を、半分の解像度のスケーリング関数とウェーブレット関数による基底展開の和に分解する。1989年に Stephane G. Mallat が発表した。
本来は異なる物だが、Mathematica や MATLAB をはじめとして、多くのソフトウェアでは多重解像度解析の事を離散ウェーブレット変換と呼んでいる。離散ウェーブレット変換の本来の定義は、離散ウェーブレット変換の項目を参照。
概要.
関数 formula_1 をスケーリング関数 formula_2 で展開した上で、
下記のウェーブレット関数 formula_4 への展開を用いて、
関数 formula_1 を異なる解像度（レベル）のウェーブレット関数に展開していく手法を多重解像度解析という。
下記関数集合が基底関数として使われる。
スケーリング関数とウェーブレット関数の間にはトゥースケール関係が成立する事が必要である。
トゥースケール関係.
以下の関係が成立する場合、トゥースケール関係と呼ぶ。
formula_11 の値はウェーブレット関数ごとに異なる。例えばハールウェーブレットの場合は、formula_12 。
トゥースケール関係が成立していると、下記の式が成立する。
formula_14 を分解数列と呼ぶ。例えば、直交ウェーブレットの場合、formula_15 。
分解数列が分かると、formula_16 という変形が可能になり、これを再帰的に繰り返すと多重解像度解析になる。
2次元.
2次元の場合は、まず、関数 formula_17 をスケーリング関数 formula_2 で基底展開する。
formula_20 に対して、
合わせて、
の4つに分解する。そして、formula_28 を同じように再帰的に分解していく。
結果として、m 回繰り返すと、下記の関数集合が基底関数となる。
3次元以上も同じ。
フィルタバンクによる表現.
1レベルの変換.
信号formula_30の離散ウェーブレット変換は、一組のフィルターを通すことによって計算される。ここではリフティングスキームではなく1989年に発表された手法を説明する。
下記の関係性を満たす formula_31 と formula_32 を formula_33 から計算するのが離散ウェーブレット変換である。j は入力依存の何らかの整数。formula_34 はスケーリング関数、formula_35 はウェーブレット関数。
最初に入力信号列から formula_33 を計算しないといけないが、ハールウェーブレットや Bior2.2 双直交ウェーブレットの場合は、入力信号の解像度と同じ解像度のスケーリング関数を使えば、formula_33 = 入力信号列になり、特に計算は不要となる。
信号は、インパルス応答がformula_39であるローパスフィルタと、同じくformula_40であるハイパスフィルタを利用して分解される。得られた出力は、ハイパスフィルタから得たものを詳細係数(detail coefficients)、ローパスフィルタから得たものを近似係数(approximation coefficients)とよぶ。これら2つのフィルタは互いに密接な関係があり、直交ミラーフィルタとして知られたものである。
次に、半分にダウンサンプリングを行う。この分解によって、時間解像度は元の信号の半分になり、周波数解像度は2倍となる。これは、ハイゼンベルクの不確定性原理を満たしている。
ダウンサンプリングのオペレータとしてformula_43を使う。
以上の式をまとめると以下のようになる。
しかしながら、formula_47の計算の後にダウンサンプリングがあるため、計算に無駄がある。リフティングスキームは、この点を改善している。
多重解像度解析.
この分解は、近似係数を再び分解することによって繰り返され、これを多重解像度解析と呼ぶ。これは、異なる時間周波数の局在性を持った部分空間を枝とする二分木によって表すことが出来る。これはフィルタバンクとして知られているものである。
各々のレベルにおいて、低周波と高周波に分解される。formula_48 の長さの入力信号は、ハールウェーブレットの場合は formula_49 のレベルに分解される。他のウェーブレットの場合も、おおよそ formula_49 になるが、多少ずれる。
ハールウェーブレットの場合での、32サンプルの信号を例にとる。周波数レンジは0から formula_51 であり、3レベルまで分解するとすると、以下の表のようになる。
リフティングスキーム.
1996年に Wim Sweldens が新しい離散ウェーブレット変換の計算方法である（）を発表した。信号を偶数番目の要素と奇数番目の要素に分け、偶数番目の要素で奇数番目の要素を予測し、予測とのズレを高周波成分とし、残差を低周波成分とする手法である。論文の中で、自由に双直交ウェーブレット関数を定義できること、計算が高速化する事を長所に挙げている。ハールウェーブレットと双直交ウェーブレットを扱える。
定義.
閉部分空間の系列 formula_52 が次の条件1〜5を満たすとき、formula_53 は多重解像度解析を成すという。formula_54はformula_55 の閉包を表す。
最後の条件はより厳しい条件として、
が用いられる事も多い。ここではこの条件5'を用いる。
条件1は空間が入れ子になっていることを意味する。条件2は formula_73 を空間 formula_74 への直交射影作用素とすると、全ての formula_75 に対して formula_76 であることを意味する。条件3は条件1の全ての空間がスケール変換で得られる事を意味する。条件4は平行移動に対して空間が普遍であることを意味する。
閉部分空間の集まりが条件1〜5'を満たすときには、いつでも正規直交ウェーブレット formula_77 による formula_78 の基底 formula_79 が存在し、formula_80 が成立する。
プログラム例.
Python（ハールウェーブレット）.
ハールウェーブレットの離散ウェーブレット変換は下記のの数式にて行える。入力を x として、formula_81 は入力 x の先頭を1番目として偶数番目の要素、つまり、codice_1、formula_82 は奇数番目の要素、つまり、codice_2。c が低周波成分（基底関数がスケーリング関数）、d が高周波成分（基底関数がウェーブレット関数）。
ハールウェーブレットの離散ウェーブレット変換のソースコードは下記のようになる。入力は x で NumPy の配列で与える。多重解像度解析をしたい場合は、codice_3 して長さが1になるまで繰り返す。formula_84 の形の基底関数の係数にしたい場合は、codice_4 と codice_5 をする。
d = x[0::2] - x[1::2]
c = x[1::2] + d / 2
ハールウェーブレットの逆離散ウェーブレット変換のソースコード。
x[1::2] = c - d / 2
x[0::2] = d + x[1::2]
Python（Bior2.2双直交ウェーブレット）.
Bior2.2双直交ウェーブレット（2階Bスプライン、線形スプライン）の離散ウェーブレット変換は下記のリフティングスキームの数式にて行える。formula_85 は一つ次の要素、formula_86 は一つ前の要素を表す。
Bior2.2の離散ウェーブレット変換のソースコード。配列の境界で足りない分は対称パディングで埋めている。
y = np.pad(x, (4, 2), "symmetric")
d = y[2::2] - (y[1:-2:2] + y[3::2]) / 2
c = y[3:-2:2] + (d[:-1] + d[1:]) / 4
d = d[1:]
Bior2.2の逆離散ウェーブレット変換のソースコード。
x[1::2] = c[1:] - (d[:-1] + d[1:]) / 4
x[0] = 2 * d[0] + x[1]
x[2::2] = d[1:-1] + (x[1:-2:2] + x[3::2]) / 2
Java.
もっとも単純な例である。ハールウェーブレットによる多重解像度解析をJavaで記述した。
/**
 * 注意：このメソッドは input 配列を破壊する。
 * input.length は2以上の2のべき乗でないといけない。
public static int[] calc(int[] input) {
 int[] output = new int[input.length];
 // length は 2^n で n は1つずつ減っていく
 for (int length = input.length » 1; ; length »= 1) {
 for (int i = 0; i < length; i++) {
 int a = input[i * 2];
 int b = input[i * 2 + 1];
 output[i] = a + b;
 output[i + length] = a - b;
 if (length == 1)
 return output;
 // 次の反復のために配列を入れ替える
 System.arraycopy(output, 0, input, 0, length « 1);
下記コードは上記を逆ウェーブレット変換する。
 * 注意：このメソッドは output 配列を破壊する。
 * output.length は2以上の2のべき乗でないといけない。
public static int[] inverse(int[] output) {
 int[] input = new int[output.length];
 for (int length = 1; ; length *= 2) {
 for (int i = 0; i < length; i++) {
 int a = output[i];
 int b = output[i + length];
 input[i * 2] = (a + b) / 2;
 input[i * 2 + 1] = (a - b) / 2;
 if (length == output.length » 1)
 return input;
 // 次の反復のために配列を入れ替える
 System.arraycopy(input, 0, output, 0, length « 1);
C言語.
C言語による、 9/7 ウェーブレット変換（JPEG-2000で利用）のリフティングを用いた高速な実装は、dwt97.c を参照。