論 文
分母分離構造を用いた逆伝達縦続法による
2次元再帰形適応フィルタ
戸島國仁 大木真 周欣欣 橋口住久
（平成7年8月31日受理）
2-D Backpropagation Cascade Adaptive Recursive Filters
with the Separable Denominator Function
KunihitoTOSHIMA MakotoOHKI XinxinZHOU SumihisaHASHIGUCHI
Abstract In this paper，1−D backpropagation cascade adaptive filters are extended to the 2−Dcase． The separable denominator function is used ill our approach． This makes it easy to perform computations for updating filter coefHcients and check the stability of adaptive filters． The system identification and the noise reduction of a image are shown to illustrate the performance of this approach．
1 まえがき
従来適応フィルタは1次元ディジタル信 号処理を対象に，零点だけが調整される非 再帰形適応フィルタや，極と零点双方が調 整される再帰形適応フィルタが研究され発 展してきた1）．最近では，1次元適応フィル タの理論を2次元信号処理の分野に拡張し た2次元適応フィルタが研究され，画像強 調や画像の雑音除去などへの応用が期待さ れている2）3）．1992年にGaoらは，逆伝達縦続法と呼
ばれる1次元再帰形適応フィルタの新しい アルゴリズムを提案した4）．彼等の方法は， 再帰形適応フィルタの伝達関数において分 母多項式を因数分解することにより再帰部 分を低次フィルタの縦続構造とし，目的信 号を再帰部分の逆フィルタに逆伝達させる 方法である．この方法により適応アルゴリ ＊電子情報工学科，Department of Electrical Engi− neering a皿d Computer Science ズムの計算量は低減され，システムの安定性判別も簡単になる．本論文ではGaoら
のアプローチを2次元再帰形適応フィルタ の場合に拡張する．2次元再帰形適応フィ ルタに拡張する場合，2通りの方法がある． 1つは2次元再帰形適応フィルタの再帰部 分を低次の2次元フィルタの縦続構造とす る方法である．この方法は棟安らにより提 案されている7）．彼等の方法では，伝達関数の分母多項式を1×1次の2次元フィル
タの縦続構造とし，目的信号をその逆フィ ルタへ逆伝達させ適応アルゴリズムを計算 している． 本論文では，もう1つの方法である分母 分離構造を用いた逆伝達縦続法を提案する． 本方法では，再帰部分を1次もしくは2次 の水平方向の1次元フィルタと1次もしく は2次の垂直方向の1次元フィルタの縦続 構造にし，目的信号を逆伝達させる．適応アルゴリズムとしてLMSアルゴリズムを
用い，再帰部分は1次元フィルタの縦続構造であるため，フィルタリングとLMSア
ルゴリズムにおける計算量を低減させるこ とができる．また，安定性判別は1次元の 安定性テストに基づき，簡単化させること ができる．シミュレーションとして，2次 元システム同定と画像の雑音除去を行い， その有効性を示す．2 逆伝達縦続法による2次元LMS
アルゴリズム分母分離構造による2次元再帰形適応
フィルタの伝達関数は次式で表わすことが できる． B（Zl，12） H（Zl，z2）＝ A・（Z・）A2（Z2） （1） 上式においてB（Zl，Z2）は21とz2の2変数 多項式であり，Al（Zl）と、42（z2）はそれぞ れ水平方向と垂直方向の1変数多項式であ る．Al（z1）とA2（z2）は1次もしくは2次の 1変数多項式に因数分解することができる． H（Zl，2「2）＝ B（Zl，212） n5＝，．41た（Zl）n品．42，（Z2） （2） り，y（m，n）は適応フィルタからの出力信 号である．d（m，n）は目的信号であり，適応 フィルタの再帰部分の逆フィルタに逆伝達 される．中間誤差信号は次式で表わされる． ei（m， n）＝di（m， n）一蝋m， n） （i＝1，… ，K−←1ン） （6） di（m， n）は目的信号を再帰部分の各逆フィ ルタに通して得られる信号であり，Yi（m，n） は適応フィルタの中間出力である． 適応フィルタの非再帰部分B（Zl，Z2）の 係数bpq（m， n）は2次元LMSアルゴリズ ムを使うと次式により更新される2）． ただし， P （？ B（・・，・・）＝ΣΣb，，・rP・；q P＝Og＝0 A、k（・、）＝1−・、ん、・r1−・、A2・r2 A21（z2）＝1一α21、・∫’一α212ぢ2 （3） （4） （5） である．直接構造による2次元再帰形フィ ルタに比べ，再帰部分は低次1次元フィル タの縦続構成であるため，フィルタリング における計算量を低減させることができる． 逆伝達法では，目的信号を逆伝達させ， 中間誤差信号を発生させる．その中間誤差 信号の2乗平均値を最小にするようにフィ ルタ係数を調整する．ブロック図を図一1に 示す．図一1において，m，nをそれぞれ水平 方向と垂直方向のインデックスとすると， u（m，n）は適応フィルタへの入力信号であ bpq（m＋1，n） ∂el（m，n） ＝bpq（m， n）一μb ∂bpq（m，n） ＝bpq（m，n）十2μbe1（m， n）u（m−P， n−9） （7）水平方向の1次元逆フィルタAlk6Dフィ
ルタ係数alks（m，n）（s＝1，2）と垂直方 向の1次元逆フィルタA216Dフィルタ係数 α21t（m，n）（t＝1，2）は1次元LMSアルゴ リズムにより更新される． αlk。（m＋1，n） ∂e2（m，n） ＝α1k，（m，n）一μ。1 ∂α1兎，（7η，n） ＝α・ん、（m，n） 十2μ。1ek（m，n）dk（m−s，n） （3＝＝1，2 k＝1，… ，K） （8） α21t（m＋1，n） ∂eえ＋1（m，n） ＝α21t（m， n）一μα2 ∂α21t（m，n） ＝α21t（m，n） 十2μα2e、K＋1（m， n）dk’＋1（m， n−t） （t＝1，2 1＝1，… ，1ン） （9） μb，μ、1，μ。2はステップサイズパラメータと 呼ばれ，収束状況を制御するパラメータでu（m，n） yκ＋・（m，n） ・41κ（Zl） eκ（m，n）i ・41κ（z・） dK（m， n） A2・（z2） κ＋・（m，n）i A2・（X2） dK＋・（m， n） dκ＋L（m， n） y（m，n） d（m，n） 図一1 逆伝達縦続法を用いた2次元分母分離再帰形適応フィルタ Fig．1 The 2−D separable denominator ADF with backpropagation methOd． ある．係数を更新することにフィルタ係数 alks，a21tは逆フィルタから適応フィルタの 再帰部分ヘコピーされ，適応フィルタの再 帰部分の係数を更新することになる． 中間誤差信号を用いることにより，出力 誤差を用いてフィルタ係数を更新するLMS アルゴリズムよりもはるかに計算を簡単に することができる．
3 安定性判別
分母分離構造による2次元再帰形適応
フィルタの利点として，簡単な安定性判別 が挙げられる．再帰部分は低次の1次元フィ ルタの縦続構造であるため，1次元フィル タの安定性チェックをそのまま用いること ができる． 再帰部分の伝達関数が式（4）で表わされ るとき，次式の条件が満たされる場合シス テムは安定であることが知られている． 1＋αlkl 一・alk2＞0 （10） 1一α・k・一α・k2＞0 （11） 1十alk2＞0 （12） 式（5）の場合も同様である．したがって図 一2に示される三角形内部にフィルタ係数が 存在するとき適応フィルタは安定である．＿1
1安 α1先2，α2ε2 @不安定領域 @ 1領域 α1丘1，α2∼1＿1
図一2 安定領域 Fig．2 The stability region．4 シミュレーション
4．1 2次元システム同定
シミュレーションとしてまず始めに2次 元システム同定を行った．システム同定の ブロック図を図一3に示す．適応フィルタは 誤差信号を減少させるようにフィルタ係数 を調整し，未知システムの伝達関数を推定 する． 本来未知である2次元システムとして，入力信号umn
L 」出力信号y（姻
2次元適応フィルタ 図一3 システム同定のブロック図 Fig．3 Block diagram of the system iden− tification． ここでは原点において円対称な4×4次の 低域通過型フィルタを用いた．伝達関数を 式（13）に示す H（・・…）一dぴ綴±）（・3）
ただし， 凡（Zl，Z2） ＝1＋・r2＋・∫2＋・「2・r2 ＋bi、（・r1＋ぢ’＋・r2・；1＋・「1万2） ＋bi、・「1・r1 Di・（z、）＝1＋・i、・r1＋・i2z「2 Di2（z2）＝1＋・i、ヨ1＋・i2zi2 A＝0．625858×10−18 （14） （15） （16） α11＝−0．216192 α12＝0．507798 α21＝−0．565076 α22＝0．183265 b11＝0．431548×109 b12＝−0．207647×109 b21＝0．131790×109 b22＝0．149429×109 である．入力信号は平均0分散1の白色ガウス
雑音とした．ステップサイズパラメータは μα1，μ。2，μbともに0．03とした．フィルタ 係数の初期値はすべて0とし，不安定な領 域に係数が更新される場合はその回の係数 更新を行なわないようにした． 一20 （−40 自 首一60 9◎ −80 一100 0 100 200 300 400 500 Number of Iterations（×1000） 図一4 平均2乗誤差の収束状況 Fig．4 Convergence curve of MSE． シミュレーション結果を図一4に示す．図 一4において横軸は係数更新の回数であり， 縦軸は目的信号と適応フィルタの出力信号 との平均2乗誤差である，係数更新回数が 約500000回で平均2乗誤差は約一100dBに なり，更に収束していく傾向にあることが 分かる． 4．2 画像の雑音除去 次に，適応ラインエンハンサによる画像 の雑音除去を行った．適応ラインエンハン サのブロック図を図一5に示す．図一5におい てx（m，n）は原画像の信号であり，ω（m，n） は雑音である．d（m，n）は目的信号であり， 原画像x（m，n）に雑音ω（m，n）を加えた信 号である．2次元再帰形適応フィルタへの入 力信号u（m，n）は目的信号を遅延させた信 号であり，過去の入力信号から現在の入力 信号を予測する予測器として適応フィルタ は動作する．目的信号を予測する際に，自 己相関の強い画像信号は予測しやすく，自 己相関の低い雑音は予測しにくいため，雑 音は低減され画像は強調される． ここでは適応フィルタへの入力信号とし て，目的信号を水平方向，垂直方向共に1 画素つつ遅延させた信号を用いた． u（m，n）＝d（m−1，n−1） （17） 原画像として標準画像データベースSIDBA のなかから大きさ256×256で256階調の mandril1を用いた．雑音は平均0分散20004mlη ω（m，n） m，η1 遅延 姻π、次元再帰形 y（m，π 十 適応フィルタ 一 e（mtn） 図一5適応ラインエンハンサのブロック図 Fig．5 2−D adaptive五ne enhancer． の白色ガウス雑音とした．2次元再帰形適 応フィルタは4×4次のものを用い，フィ ルタ係数の初期値は全て0とした．また， 出力画像の評価を行なうためp．pSNRを次 式により計算した． P−pSNR＝ 1010910 2552×2562 255 255 ΣΣ｛x（m，n）−y（m，n）｝2 m＝On ＝O （18） 図一6に原画像，図一7に雑音を加えた画像を
示すこのときのp．pSNRは15．1dBであ
る．図一8は文献7）により提案された2次 元再帰形適応フィルタによる出力画像であ り，伝達関数の分母多項式を低次の2次元 多項式の縦続構造にして目的信号を逆伝達 させる方法である．図一9は本研究で提案し た2次元再帰形適応フィルタによる出力画 像である．図一7において白い斑点のように 見える雑音が図一8，9において除去されて いるのが分かる． 表1は文献7）で提案された方法と本研究 による方法とのステップサイズパラメータとp．pSNRである．表1からp−pSNRは2
つの方法ともほぼ等しいことが分かる．表 2は本研究による方法と文献7）で提案され ている方法との1回のフィルタリングに要 する計算量の比較である．本研究では分母 分離構造を用いているためフィルタリング に要する計算量が文献7）のものよりも少 なく，しかもほぼ同等の雑音除去能力があ ることが分かる．5 まとめ
本論文では分母分離構造を用いた逆伝達縦続法による2次元再帰形適応フィルタ
を提案した．この方法ではフィルタリング に要する計算量を低減させることができ， フィルタ係数の更新式と安定性判別も簡単 にすることができることを示した． シミュレーションとして2次元システム 同定と画像の雑音除去を行い，提案した方 法が妥当であることを示した．また画像の 雑音除去では，分母多項式を2次元の低次 多項式の縦続形にしたものとの比較を行い， ほぼ同等な結果を得られることを示した．参考文献
1）J．J．Shynk： ‘‘Adaptive IIR Filter− ing，，， IEEE ASSP Magazine． pp．4− 21，April，1989． 2）M．M．Hadhoud and D．W．Thomas： “The Two−Dimensional AdaptiveLMS（TDLMS）Algorithm”， IEEE
Trans． Circuits Syst．， Vol．35， No．5， pp．485−494， May，1988． 3）M．Ohki and S．Hashiguchi：“Two− dimensional LMS adaptive filters，，， IEEE Trans． Consum． Electron．， Vol．37， No．1， pp66−73， Feb，1991． 4）F．X．Y．Gao and W．M．SnelgroVe：“An Adaptive Backpropagation Cascade IIR Filter，，， IEEE Trans Circuits． Syst．，39，NO．9，pp．606−610，Sep，1992． 5）吉江智隆，川又政征，樋口龍雄：“並 列完全分離形2次元適応フィルタの実 現”，第8回ディジタル信号処理シン ポジウム講演論文集，pp．159−165，1993 年10月． 6）棟安実治，雛元孝夫：“LMS法を用 いた2次元適応状態空間フィルタ”， 電子情報通信学会論文誌，Vol．J．76− A，No．4，pp．598−604，1992年4月． 7）棟安実治，上本栄治，石崎直哉，雛元 孝夫：“逆伝搬法を用いた2次元縦続形 IIRフィルタ”，1994年電子情報通信 学会春季大会予稿集，A−167，1994年図一6 原画像 Fig．6 0riginaユimage． 図一7 雑音を加えた画像 Fig．7 1mage degraded l）y Gaussian noise． （p−pSNR＝15．1 dB） 図一8 文献7）の方法による出力画像 Fig．8 0utput image by Muneyasu’s method． （p−pSNR＝18．3dB） 図一9 提案した方法による出力画像 Fig．9 0utput image by the proposed method． （p．pSNR＝18．1 dB） 表一1 出力画像の比較 Table−1 p−pSNR of resultant images． 文献7）の方法 提案した方法 μ 3．0×10−8 3．0×10−8
