64-4-1
サンプル値 H ∞ 制御理論による最適ウェーブレット展開○小蔵正輝 永原正章 山本裕 （京都大学）

Optimal Wavelet Expansion via Sampled-Data H ∞ Control Theory∗M. Ogura, M. Nagahara, and Y. Yamamoto (Kyoto University)
Abstract– Wavelet expansion of an L2 signal requires the L2 inner product of the original signal and ascaling function. This process is non-causal because the integral is to be executed over the whole time interval(−∞, ∞). In digital signal processing, it is common to use sampled values of continuous-time signals insteadof the inner product. This however follows a large reconstruction error, called ”wavelet crime.” We therefore
design a causal system which produces an approximation of the inner product via sampled-data H ∞ control
theory. We treat not only the Haar wavelet but also general wavelets including the Daubechies wavelet when
the scaling function has a finite support. We also consider a time-varying model of input signals.
Key Words: Wavelet expansion, Sampled-data control, Signal reconstruction.

1

はじめに

連続時間信号を時間と周波数の両面から解析する手法として，ウェーブレット展開がある．これは，時間と周波数に局在した基底で信号を展開するという手法
であり，その展開係数から信号の局所的な振る舞いを読み取ることができる．
連続時間信号に対してウェーブレット展開を行うには，信号とスケーリング関数の内積を計算する必要があり，信号の任意の時刻における値が必要である．しかし，連続時間信号の全ての値がわかることは現実にはなく，一般にはそのサンプル値しかわからないことがほとんどである．この場合，ウェーブレット展開を直接行うことはできない．
これと全く同じ問題が，ウェーブレット変換についても存在する．そのため，サンプル値からウェーブレッ
ト係数を計算するための手法がいくつか提案されてきた．最も簡単なものとして，信号のサンプル値をそのままウェーブレット係数として用いるという便宜的な手法がある．これは単純さと引き替えに大きな誤差を生じさせるため，“wavelet crime” と呼ばれている 8) ．
サンプル点間を標本化関数 sinc := sin(π •)/(π •) で補間することにより内積を計算するという手法 1) や最適化による手法 3, 5, 7) なども提案されているが，いずれの手法も未来における信号のサンプル値を全て必要とするため，非現実的である．有限長のフィルタによるサンプル値からのウェーブレット係数の計算手法としては 4) があり，サンプル値制御理論を用いて，信号からウェーブレット係数の計算誤差への作用素の誘導ノルムを最小化するフィルタを設計している．
ウェーブレット展開に戻ると，この展開はある解像度で実際には打ち切られる．展開する信号の解像度が，この打ち切られた解像度に制限されている場合には，ウェーブレット係数の計算誤差はウェーブレット展開の誤差にそのまま反映されるため，先に挙げた既存の手法のようにウェーブレット係数を精度良く計算することを考えればよい．しかし，信号の解像度が制限されていることは，現実にはあり得ない．この場合，展開誤差を小さくするためには，信号の高解像度成分についても考慮して展開係数を定める必要がある．
そこで本研究では，サンプル点間を含めて展開誤差

第7回制御部門大会（2007年3月6日〜8日・調布）を最小にするような最適ウェーブレット展開係数を信号のサンプル値から計算する手法を提案する．サンプル
値 H ∞ 制御理論の枠組みで問題の定式化を行い，H ∞最適な離散時間コントローラの設計に帰着することを示す．そして，数値計算例を用いて既存手法に対する提案法の優位性を示す．また，信号の時変モデルに対する適用についても考察する．

2

ウェーブレット理論

2.1 多重解像度解析とスケーリング関数
まず，ウェーブレット理論において重要な役割を果たす多重解像度解析とスケーリング関数を定義する 6) ．
定義 1 (多重解像度解析). 空間 L2 (R) における閉部分
空間列 {Vj }j∈Z で，次の条件 1〜5 を満たすものを，多重解像度解析という．さらに，条件 5 に現れる φ をス
ケーリング関数という．

1. ∀j ∈ Z, Vj ⊂ Vj+1 .
TS
2. j∈Z Vj = {0}, j∈Z Vj = L2 (R).
3. ∀j ∈ Z, f ∈ Vj ⇔ f (2 •) ∈ Vj+1 .

4. ∀k ∈ Z, f ∈ V0 ⇔ f (• − k) ∈ V0 .

5. ある φ ∈ L2 (R) が存在し，{φ(• − k)}k∈Z が V0 のRiesz 基底である．
今，f ∈ Vj の Fourier 変換の台 Ω := supp fˆ が有限であるとすると，f (2 •) ∈ Vj+1 の Fourier 変換の台は Ωを 2 倍に拡げた 2Ω となる．よって，条件 3 より，Vj+1
は Vj より一段階高い周波数成分を含むと解釈できる．
この高周波成分を生成する差分である，Vj の Vj+1 における直交補空間を Wj := Vj+1 ª Vj と定義する．
すると条件 2 は，j の増加と共に Vj に属する信号が細かく，即ち解像度が高くなっていき，j → ∞ での極限で Vj は全空間 L2 (R) となること，また，逆に j が減少していくと Vj に属する信号の解像度が低くなって
いき，j → −∞ の極限においては恒等的に 0 の関数しか表現できなくなることを意味していると解釈できる．
スケーリング関数は各部分空間 Vj , Wj の基底を生成する．まず，平行移動 {ψ(• − k)}k∈Z が W0 の Riesz 基

SY0003/07/0000-6441 (c) 2007 SICE

底となるような ψ ∈ L2 (R) を φ から構成することができ，この ψ をウェーブレットという．さらにスケーリング関数とウェーブレットの伸縮・平行移動を

cj

q

φj,k := 2j/2 φ(2j • − k), ψj,k := 2j/2 ψ(2j • − k)で定義する．条件 3 と 4 より，{φj,k }k∈Z と {ψj,k }k∈Zがそれぞれ Vj と Wj の Riesz 基底であるとわかる．
以上より，空間 L2 (R) の直和分解が可能である:

∞M
L2 (R) = Vj0 ⊕ 
Wj  .

f=X

cj0 (k)φj0 ,k +∞ XX

dj (k)ψj,k

2.2 ピラミッドアルゴリズム
ある解像度のウェーブレット係数からそれより低い
解像度のウェーブレット係数を高速に求める方法にピ
ラミッドアルゴリズムがあり，これにより信号の高速な解析が可能である．簡単のため，スケーリング関数の平行移動が V0 の正規直交底を構成する場合を考える．このとき，{φj,k }k∈Z と {ψj,k }k∈Z はそれぞれ Vjと Wj の正規直交底となるので，f ∈ L2 (R) のウェー
ブレット係数は次のように内積で与えられる．
cj (k) = (f, φj,k ) , dj (k) = (f, ψj,k ) .

(1)

また，φ ∈ V0 と ψ ∈ W0 を V1 の基底関数でそれぞれ
次のように展開することができる．
XX
φ=hn φ1,n , ψ =gn φ1,n .
(2)
n∈Z

n∈Z

式 (1)，(2) より，解像度 j の係数から解像度 j − 1 の係数を次のように計算することができることがわかる．
X
cj−1 (k) =hn−2k cj (n),n∈Z

dj−1 (k) =X

gn−2k cj (n).

(3)

n∈Z

さらに，この計算は Fig. 1 のマルチレートフィルタバ
ンクにより高速に行われる．ただし，H と G はインパルス応答がそれぞれ {h−n }, {g−n } であるフィルタ，↓ 2
はダウンサンプラ (↓ 2y)(k) := y(2k) である．Fig. 1 のフィルタバンクの直列結合を用いれば，ある解像度のウェーブレット係数からそれより低い解像度のウェー
ブレット係数を高速に計算することができる．この計算法をピラミッドアルゴリズムという．
ピラミッドアルゴリズムを用いると，解像度制限されている関数のウェーブレット展開は，その関数の持つ最高解像度のウェーブレット係数を求めることに帰

- H

- ↓2

cj−1
-

X

cj0 (k)φj0 ,k +J−1XX

dj (k)ψj,k

j=j0 k∈Zk∈Zとなり，この展開に現れる各係数は，係数 cJ のみからピラミッドアルゴリズムにより計算される．

3と各部分空間の基底で展開できる．これをウェーブレッ
ト展開といい，信号 f から係数 {cj0 , dj0 , dj0 +1 , . . . } への写像を離散ウェーブレット変換という．

dj−1

着される．今，f ∈ VJ とすると，直和分解 VJ = Vj0 ⊕LJ−1( j=j0 Wj ) に対応して f のウェーブレット展開はf=

j=j0 k∈Z

k∈Z

- ↓2

Fig. 1: 式 (3) に対応するマルチレートフィルタバンク

j=j0

さらに，この直和分解に対応して任意の f ∈ L2 (R) を- G

最適化問題への帰着

この章では，最適なウェーブレット展開係数を信号のサンプル値から計算する離散時間コントローラの設計問題を定式化し，これが離散時間の H ∞ 最適化問題に帰着されることを示す．

3.1 設計問題の定式化
第 2 節で述べたように，信号の解像度制限仮定のもとで，ウェーブレット展開はその信号の持つ最高解像度のウェーブレット係数を求めることに帰着される．しかし，周波数帯域制限された信号が現実には存在しないのと同様に，解像度が制限された信号も現実には存在しない．また，本稿のはじめに述べたように，連続時間
信号の全ての値が利用できることも現実にはなく，そのサンプル値しか利用できないことがほとんどである．
そこで，信号のサンプル値しか利用できないという条件のもとで，サンプル点間を含めて誤差が小さい最適ウェーブレット展開係数を信号のサンプル値から計算することを考える．
問題の定式化の前に，いくつかの仮定を行う．まず，対象とするウェーブレットのスケーリング関数が有限の台を持つと仮定する．
仮定 1. ある m > 0 が存在し，supp φ ⊂ [0, m].

(4)

また，展開の対象とする信号のクラスを L2 [0, ∞) とし，信号のサンプリング周期は 1 に正規化する．さらに，展開を行う解像度は 0 とする．これらの仮定によって一般性は失われない．
そして，Fig. 2 の誤差系を考える．ただし，W はそ

- Dl
w -

W

q- S

- Kd

c-

+?
de H̃φ −Fig. 2: 誤差系の伝達関数が厳密にプロパーかつ安定である連続時間
線形時不変システムであり，対象とする信号の周波数特性を表す．Dl は整数 l だけの連続時間における遅れ
Dl f := f (• − l) であり，(4) を満たす最小の m より大き

くとることにより展開の因果性を確保する必要がある．
Kd は設計する離散時間コントローラである．最後に，H̃φ は逆離散ウェーブレット変換であり，解像度 0 における，入力の離散時間信号を展開係数としたウェーブ
レット展開を出力する．

H̃φ : `2R → L2 [0, ∞) : {c(k)}∞k=0 7→∞
Xc(k)φ0,k .

k=0

φ がリース基底であるため，H̃φ は確かに `2R からL2 [0, ∞) の中への写像となっている．そして，次のよ
うに問題を定式化する．
問題 1. Fig. 2 の誤差系における w から e への作用素のL2 [0, ∞) 誘導ノルムを最小にする離散時間コントロー
ラ Kd を求めよ．

Fig. 3: 誤差系の周波数応答：実線 提案法, 鎖線 wavelet
crimeFig. 2 の 誤 差 系 を 下 側 線 形 分 数 変 換 を 用 い てFl (G, Kd ) と表すことができる．ただし，·¸
Dl W −H̃φ
G :=
.
SW
0
すると，問題 1 は誤差系の高速サンプル／ホールド近似

Fl (Gd , Kd ) := LN S1/N Fl (G, Kd )H1/N L−1
N(5)の `2RN 誘導ノルム最小化問題で近似できる 10) ．ただし，
S1/N は周期 1/N の理想サンプラ，H1/N は周期 1/Nの 0 次ホールド，LN は N 次の離散時間リフティングである．そして，問題 1 の近似として次の問題を解く．

Fig. 4: 復元信号の時間応答：実線 提案法，鎖線 wavelete
crime，点線 原信号，白丸 サンプル点問題 2. kFl (Gd , Kd )k を最小にする，1 入力 1 出力離散時間コントローラ Kd を求めよ．
式 (5) の Gd は有限次元の離散時間線形時不変システ
ムであるから，この問題は離散時間の H ∞ 最適化問題
であり，MATLAB 等を用いて解くことができる．

4

数値例

この章では，数値計算例を示すことにより提案法の有効性を示す．対象とする信号の周波数特性 W の伝達関数 W は次のように定めた．

W (s) :=

1.
(10s + 1)(s/3 + 1)

まず提案法と wavelet crime を比較する．スケーリ
ング関数を有限の台を持つ Daubechies の直交ウェーブ
レット 2) のスケーリング関数（バニッシングモーメン
トは 2 とする），許容する遅れ l を 3，高速サンプル／
ホールドの次数 N を 8 としてコントローラを設計した．
Fig. 2 の誤差系におけるコントローラ Kd に，提案法によるコントローラと wavelet crime によるコント
ローラ（1 の遅延）をそれぞれ代入した誤差系のサン
プル値制御系としての周波数応答 9) が Fig. 3 である．
鎖線の wavelet crime と比べ，実線の提案法は全域的に誤差が小さい．
また，三角波信号のサンプル値からの復元を両手法により行った．Fig. 4 は復元信号の時間応答，Fig. 5 は復元誤差の時間応答であり，実線の提案法による誤差
は鎖線の wavelet crime による誤差より小さい．

Fig. 5: 復元誤差の時間応答：実線 提案法，鎖線 waveletcrime
次に，加嶋らによる従来法 4) との比較を行う．スケーリング関数は 4 階の B-spline とし，許容する遅れ l を
6，高速サンプル／ホールドの次数 N を 8 としてコントローラを設計した．また，比較の対象して，従来法を用いて 19 の遅れをもつコントローラを設計した．
それぞれのコントローラを代入した誤差系の周波数応答が Fig. 6 であり，実線が提案法，破線が従来法である．提案法は許容する遅れが半分以下であるにもかかわらず，従来法より良い性能を示している．
また，矩形波信号のサンプル値からの復元を両手法により行った．Fig. 7 は復元信号の時間応答，Fig. 8 は復元誤差の時間応答である．従来法による復元信号に生じているリップルの大きさが，提案法では小さくなっている．また，誤差も小さくなっている．

5

信号の時変モデル

3 節の定式化では，信号のモデルとして線形時不変システムによる周波数特性仮定を用いた．この拡張と

今，Ḡd が有限次元であることは簡単に確かめることができる．故に，この問題は W が時不変である場合の問題 2 と同様，離散時間の H ∞ 最適化問題である．得
られるコントローラ K̄d の入出力は次元 T を持つが，この入出力に T 次の逆離散時間リフティングを施して更に遅延をかけることにより，周期 T を持つ因果的なT −1
最適コントローラ Kd := L−1
を問題 3 のT K̄d LT U
近似解として得る．ただし，U は離散時間における 1の遅延である．

6
Fig. 6: 誤差系の周波数応答：実線 提案法, 破線 従来法

終わりに

本稿では，スケーリング関数が有限の台をもつウェーブレットについて，信号のサンプル値からサンプル点間を考慮した最適なウェーブレット展開係数を計算するコントローラの設計手法を提案した．そして，コン
トローラの設計が H ∞ 最適な離散時間コントローラを設計する問題に帰着することを示し，数値計算例により提案法の有効性を確認した．最後に，信号の時変モデルに対する適用についても考察した．

参考文献

Fig. 7: 復元信号の時間応答：実線 提案法，破線 従来
法，点線 原信号，白丸 サンプル点

1) P. Abry and P. Flandrin. On the initialilzation of
the discrete wavelet transform alogrithm. IEEE
Signal Proc. Lett., Vol. 1, No. 2, pp. 32–34, Feb.
1994.
2) I. Daubechies. Ten Lectures on Waveletes. SIAM,
1992.
3) T. G. Dvorkind, H. Kirshner, Y. C. Eldar, and
M. Porat. Minimax approximation of representation coefficients from generalized samples. Technion, CCIT Report, No. 582, Apr. 2006.
4) K. Kashima, Y. Yamamoto, and M. Nagahara.
Optimal wavelet expansion via sampled-data control theory. IEEE Signal Proc. Lett., Vol. 11,
No. 2, pp. 79–82, Feb. 2004.

Fig. 8: 復元誤差の時間応答：実線 提案法，破線 従来法
して，整数周期を持つ周期時変システムによる信号のモデル化が考えられる．ここで連続時間システムが周期 T を持つことを，システムの状態方程式の各係数が周期 T を持つことで定義する．そして，問題 1 で W が
周期 T を持つ場合として次の問題を考える．
問題 3. Fig. 2 において，W が周期 T を持つとする．
このとき，w から e への作用素の L2 [0, ∞) 誘導ノルムを最小にする離散時間コントローラ Kd を求めよ．
この問題も，高速サンプル／ホールド近似により離散
時間の H ∞ 最適化問題に帰着できる．W が周期 T を持つため，式 (5) の高速サンプル／ホールド近似において
Gd は離散時間の周期 T を持つ．よって，Gd の入出力に
T 次の離散時間リフティングを施した Ḡd := LT Gd L−1
Tは時不変である．さらに，K̄d := LT Kd L−1
T とすると，Fl (Gd , Kd ) = L−1
T Fl (Ḡd , K̄d )LT が成り立つ．よって，問題 3 の近似として次の問題を得る．
問題 4. kFl (Ḡd , K̄d )k を最小にする，T 入力 T 出力離散時間コントローラ K̄d を求めよ．

5) H. Kirshner and M. Porat. On sampling invariant signal represantation. Technion, CCIT report
520, Feb. 2005.
6) S. G. Mallat. A theory of multiresolution signal
decomposition: the wavelet representation. IEEE
Trans. Patt. Anal. Machine Intell., Vol. 11, pp.
674–693, 1989.
7) P. Qian and B. A. Francis. Optimal initialization of the discrete wavelet transform. In Proc.
Workshop Recent Advances in Control, 1998.
8) G. Strang and T. Nguyen. Wavelet and Filter
Banks. Wellesley-Cambridge Press, 1997.
9) Y. Yamamoto and P. P. Khargonekar. Frequency
response of sampled-data systems. IEEE Trans.
Automat. Contr., Vol. 41, pp. 166–176, 1996.
10) Y. Yamamoto, A. G. Madievski, and B. D. O.
Anderson. Approximation of frequency response
for sampled-data control systems. Automatica,
Vol. 35, No. 4, pp. 729–734, 1999.