JP2016219860A - 未知伝達系推定装置、未知伝達系推定方法、およびプログラム - Google Patents

Abstract

【課題】適切に未知伝達系の特性を推定する。【解決手段】フィルタ係数記憶部２０１は、未知伝達系の伝達特性を表すフィルタ係数を記憶する。入力信号絶対値算出部３０１は、入力信号の絶対値の系列を算出する。出力信号推定部２０２は、入力信号の系列とフィルタ係数との畳み込み演算により推定出力信号を算出する。残差信号算出部２０３は、観測出力信号と推定出力信号との残差信号を算出する。残差信号絶対値算出部３０２は、残差信号の絶対値を算出する。フィルタ期待値算出部３０３は、フィルタ係数および外乱信号の振幅値が優ガウス性分布に従うものとして観測出力信号が観測された条件のもとでフィルタ係数が取り得る期待値を算出し、フィルタ係数を更新する。【選択図】図２

Description

この発明は、電気信号や音、振動が伝播する伝達系の特性を推定する技術に関する。

電気信号や音、振動が伝播する未知の伝達系の特性を推定することは、当該伝達系への任意の入力に対する当該伝達系の出力を予測する上で有益である。伝達系の特性を推定する問題は、有限のパラメータで表現可能なモデルをあらかじめ与え、当該パラメータの値を推定する問題として扱うことが多い。あらかじめ与えるモデルは、推定対象の伝達系の物理的素性（入力と出力の関係）に基づき選定される。この発明で対象とするのは、未知の伝達系への入力とそれに対する出力とが線形な関係を有し、その伝達系の特性がFIR（Finite Impulse Response）フィルタによってモデル化できる場合における、当該フィルタ係数の推定の問題とする。

伝達系への入力信号x(n)に対し、観測される出力信号y(n)は、式（１）で表すことができる。

ここで、h₀, h₁, …, h_Lは推定すべき未知の有限長インパルス応答の各要素である。ただし、Lは0以上の整数であり、nは離散時間を表す。また、v(n)は出力信号の観測時に混入し得る外乱信号である。

フィルタ係数W(n)は、時刻nにおけるh₀, h₁,…, h_Lの推定値を要素としてもつベクトルであり、式（２）で表される。

すると、時刻nにおいて得られる残差信号r(n)は、式（３）で表される。

このとき、フィルタ係数の推定値を式（４）により更新できることが知られている（例えば、非特許文献１参照）。

ただし、

であり、μ, εは収束特性の調整のために用いられる定数である。μはステップサイズ、εは正則化係数などと呼ばれる。また、フィルタ係数W(n)の各要素について書き下すと、式（６）となる。

式（４）の更新式の妥当性については、以下のように解釈することができる。まず、出力信号y(n)が観測された条件のもとでの、更新されるべきフィルタ係数の推定値W(n+1)についての事後確率は、ベイズの定理により、式（７）が成り立つ。

なお、式（７）では簡単のため時刻nについての標記を省略した。ここで、p(y|W)は、未知伝達系の特性がWであった場合における出力信号yが取り得る値の確率（Wについての尤度)を表し、p(W)は更新前のWが取り得る値の事前確率である。p(y|W)は、出力信号yのうち、Wに起因する未知伝達系の出力が与えられたものとして、それ以外の信号成分である外乱信号v(n)の取り得る値の確率に相当する。

外乱信号v(n)およびフィルタ係数W(n)がいずれもガウス分布に従うとすると、式（７）は、式（８）で表すことができる。

ここで、α, βは外乱信号v(n)およびフィルタ係数W(n)それぞれの分散の逆数（精度）であり、_Wは更新前の推定値を表すものとする（_・は数式中では・に下線が付された記号を表す）。また、||・||₂はベクトルのL2ノルムを表す。式（８）について、両辺の対数を取ると、

となり、Jの最大値を与えるWは式（８）の事後確率の最大値を与える。そこで、

を満足するWを計算すると

が得られる。式（４）と式（11）を対比すると、式（４）において、μ=1、ε=β/αとすれば、両者は一致することがわかる。

上述の従来技術に基づく未知伝達系推定装置２００は、図１に示すように、フィルタ係数記憶部２０１、出力信号推定部２０２、残差信号算出部２０３、入力信号二乗値算出部２０４、およびフィルタ期待値算出部２０５を有する。フィルタ係数記憶部２０１は、推定されるフィルタ係数W(n)を保持する。出力信号推定部２０２は、入力信号x(n)の系列X(n)とフィルタ係数W(n)との畳み込み演算により（外乱信号v(n)が含まれない）未知伝達系の推定出力信号W(n)^TX(n)を算出する。残差信号算出部２０３は、（外乱信号v(n)が含まれ得る）未知伝達系の観測出力信号y(n)と推定出力信号W(n)^TX(n)の残差信号r(n)を式（３）により計算する。入力信号二乗値算出部２０４は、入力信号の二乗値x(n)²の系列を計算する。フィルタ期待値算出部２０５は、入力信号x(n)の系列と入力信号の二乗値x(n)²の系列と残差信号r(n)とを用いて、フィルタ係数W(n)および外乱信号v(n)の振幅値がガウス分布に従うものとして、観測出力信号y(n)が観測された条件のもとでフィルタ係数が取り得る期待値W(n+1)を式（４）により算出し、フィルタ係数記憶部２０１に保持されたフィルタ係数W(n)を更新する。

N. J. Bershad, "Behavior of the ε-normalized LMS algorithm with Gaussian inputs", IEEE Transactions on Acoustics, Speech and Signal Processing, Vol. 35, No. 5, pp. 636-644, 1987.

式（４）に基づくフィルタ係数の推定値W(n)の更新は、外乱信号v(n)およびフィルタ係数W(n)がガウス分布に従う確率変数とみなすことで得られた式（11）の更新式と相似している。すなわち、式（４）の更新式は、外乱信号v(n)およびフィルタ係数W(n)がガウス分布に従う場合には妥当であるが、そうでない場合にはより適切な更新式が存在し得るといえる。例えば、推定すべき伝達特性が音響伝達系のインパルス応答である場合、指数的な減衰特性を有する場合が多い。この場合、各タップ係数の振幅値の分布は、零近傍への集中度の高い優ガウス性を有することが想定される。また、外乱信号についても突発的あるいは断続的な雑音や音声信号のような非定常信号である場合について、ガウス分布の仮定は必ずしも適切とはいえない。

この発明の目的は、未知伝達系のインパルス応答の振幅値や観測時に混入する外乱信号の振幅値がガウス分布の仮定に沿わない場合に、より適切に未知伝達系の特性を推定することができる未知伝達系推定技術を提供することである。

上記の課題を解決するために、この発明の未知伝達系推定装置は、未知伝達系の伝達特性を表すフィルタ係数を記憶するフィルタ係数記憶部と、未知伝達系への入力信号の絶対値の系列を算出する入力信号絶対値算出部と、入力信号の系列とフィルタ係数との畳み込み演算により推定出力信号を算出する出力信号推定部と、入力信号を入力とする未知伝達系の出力を観測することで得られる観測出力信号と推定出力信号との残差信号を算出する残差信号算出部と、残差信号の絶対値を算出する残差信号絶対値算出部と、入力信号の系列と入力信号の絶対値の系列と残差信号と残差信号の絶対値とを用いて、フィルタ係数および観測出力信号に含まれうる外乱信号の振幅値が優ガウス性分布に従うものとして上記観測出力信号が観測された条件のもとでフィルタ係数が取り得る期待値を算出し、フィルタ係数を更新するフィルタ期待値算出部と、を含む。

この発明の未知伝達系推定技術によれば、未知伝達系のインパルス応答の振幅値や観測時に混入する外乱信号の振幅値がガウス分布の仮定に沿わない場合に、より適切に伝達系の特性を推定することができる。

図１は、従来の未知伝達系推定装置の機能構成を例示する図である。 図２は、実施形態の未知伝達系推定装置の機能構成を例示する図である。 図３は、実施形態の未知伝達系推定方法の処理フローを例示する図である。

実施形態の説明に先立ち、この発明の基本的な考え方を説明する。

この発明では、従来の未知伝達系推定技術において、式（７）で定義された事後確率に優ガウス性を持つ確率分布の１つであるラプラス分布を適用し、式（12）で表される事後確率に基づいてフィルタ係数の更新値を得る。

ここで、||・||₁はL1ノルムであり、・^Tは行列の転置であり、_Wは更新前のフィルタ係数の値を要素として持つベクトルであり、αは外乱信号の精度であり、Bは非負の要素β₀, β₁, …, β_Lを対角成分にもつL+1行L+1列の対角行列である。行列Bは式（13）で表すことができる。

行列Bの対角成分が等しくβの値を取るとき、式（12）は式（14）で表すこともできる。

式（12）または式（14）で与えられる事後確率の分布は、その最大値を中心に対称な形状とならないことから、その最大値を与えるWは、期待値E(W|y)と一致しない。そこで、Wの各要素w_k（k=0, 1, …, L）について、期待値を直接計算し、式（15）として更新値を与えるものとする。

具体的には、式（12）と式（15）から、

と与える。ただし、

である。ここで、Qは、式（16）の分母と分子で相殺されうる任意の乗数であり、必ずしも必要ではない。

各フィルタ係数w_k(n)の更新値は、基本的に式（16）に基づき計算されるが、このとき用いるD_v(n), D_i(n), G_v,k(n), G_i,k(n)の計算方法は、必ずしも、式（17）から式（20）の形式と一致させる必要はなく、これらと等価あるいは近似的に等価な式に基づいて計算してもよい。あるいは、零除算を防止するなどの安定化項を加えたり、分母が小さな値を取る場合等において場合分けにより別の計算式や定数を与えたりしてもよい。さらに、式（16）の形式に限定されるものでもなく、これと等価あるいは近似的に等価な式に基づいて計算した場合についても、この発明の技術に含まれる。

以下、この発明の実施の形態について詳細に説明する。なお、図面中において同じ機能を有する構成部には同じ番号を付し、重複説明を省略する。

［第一実施形態］
第一実施形態の未知伝達系推定装置３００は、図２に例示するように、フィルタ係数記憶部２０１、出力信号推定部２０２、残差信号算出部２０３、入力信号絶対値算出部３０１、残差信号絶対値算出部３０２、およびフィルタ期待値算出部３０３を含む。この未知伝達系推定装置３００が、図３に例示する各ステップの処理を行うことにより第一実施形態の未知伝達系推定方法が実現される。

未知伝達系推定装置３００は、例えば、中央演算処理装置（CPU: Central Processing Unit）、主記憶装置（RAM: Random Access Memory）などを有する公知又は専用のコンピュータに特別なプログラムが読み込まれて構成された特別な装置である。未知伝達系推定装置３００は、例えば、中央演算処理装置の制御のもとで各処理を実行する。未知伝達系推定装置３００に入力されたデータや各処理で得られたデータは、例えば、主記憶装置に格納され、主記憶装置に格納されたデータは必要に応じて中央演算処理装置へ読み出されて他の処理に利用される。未知伝達系推定装置３００の各処理部は、少なくとも一部が集積回路等のハードウェアによって構成されていてもよい。

未知伝達系推定装置３００の備えるフィルタ係数記憶部２０１は、例えば、RAM（Random Access Memory）などの主記憶装置、ハードディスクや光ディスクもしくはフラッシュメモリ（Flash Memory）のような半導体メモリ素子により構成される補助記憶装置、またはリレーショナルデータベースやキーバリューストアなどのミドルウェアにより構成することができる。

未知伝達系推定装置３００のフィルタ係数記憶部２０１には、更新前のフィルタ係数W(n)が記憶されている。フィルタ係数W(n)は従来技術と同様であり、時刻nにおける未知伝達系１００の伝達特性h₀, h₁,…, h_Lの推定値を要素としてもつベクトルである。

図２に例示するように、未知伝達系推定装置３００および未知伝達系１００へ入力信号x(n)が入力され、未知伝達系１００は入力信号x(n)に対する出力信号を出力する。その出力信号を観測した観測出力信号y(n)が未知伝達系推定装置３００へ入力される。観測出力信号y(n)には観測時に外乱信号v(n)が混入している。未知伝達系推定装置３００では、入力信号x(n)は入力信号絶対値算出部３０１、出力信号推定部２０２、およびフィルタ期待値算出部３０３へ送られる。観測出力信号y(n)は残差信号算出部２０３へ入力される。

図３を参照して、実施形態の未知伝達系推定方法の処理手続きを説明する。

ステップＳ３０１において、入力信号絶対値算出部３０１は、入力信号x(n)の系列から入力信号の絶対値|x(n)|の系列を算出する。入力信号の絶対値|x(n)|の系列はフィルタ期待値算出部３０３へ送られる。

ステップＳ２０２において、出力信号推定部２０２は、フィルタ係数記憶部に記憶されたフィルタ係数W(n)を読み出し、入力信号x(n)の系列X(n)とフィルタ係数W(n)との畳み込み演算により推定出力信号W(n)^TX(n)を算出する。推定出力信号W(n)^TX(n)は残差信号算出部２０３へ送られる。

ステップＳ２０３において、残差信号算出部２０３は、外乱信号v(n)が含まれ得る観測出力信号y(n)から推定出力信号W(n)^TX(n)を減算した残差信号r(n)を式（３）により算出する。残差信号r(n)は残差信号絶対値算出部３０２およびフィルタ期待値算出部３０３へ送られる。

ステップＳ３０２において、残差信号絶対値算出部３０２は、残差信号r(n)から残差信号の絶対値|r(n)|を算出する。残差信号の絶対値|r(n)|はフィルタ期待値算出部３０３へ送られる。

ステップＳ３０３において、フィルタ期待値算出部３０３は、入力信号x(n)の系列と入力信号の絶対値|x(n)|の系列と残差信号r(n)と残差信号の絶対値|r(n)|とを用いて、フィルタ係数W(n)および外乱信号v(n)の振幅値が優ガウス性分布に従うものとして、観測出力信号y(n)が観測された条件のもとでフィルタ係数が取り得る期待値W(n+1)を算出する。具体的には、フィルタ期待値算出部３０３は、式（12）または式（14）で定義される事後確率を最大とするフィルタ係数の期待値を算出する。フィルタ係数の各要素の期待値は、式（16）に基づいて計算することができる。算出したフィルタ係数の期待値W(n+1)はフィルタ係数の更新値としてフィルタ係数記憶部２０１へ記憶される。

［第二実施形態］
第一実施形態では、外乱信号v(n)が優ガウス性を有するとして、式（12）または式（14）を与えた。しかし、実環境においては、ガウス性を有する定常性の強い雑音が混入することも多い。第二実施形態では、音響エコーキャンセラのように、入力信号が主に優ガウス性の強い音声信号である場合を想定することで、外乱信号がガウス性雑音であった場合にも対応できる形態を示す。具体的には、式（12）は式（21）のように修正し、式（14）は式（22）のように修正する。

ここで、f(X)は、例えば、XのL1ノルムやL2ノルムのように、入力信号の大きさに依存して非負の値を取る関数である。f(X)は、必ずしも入力信号Xのすべての要素の値に基づいて計算される必要はない。また、時刻nよりも過去の入力信号Xの値を用いたり、現在と過去の値を合わせて用いて値を算出したりすることも可能である。

この修正を実施するためには、式（16）の更新値を得るにあたり、式（17）から式（20）の計算において、αをα/f(X)で置き換えればよい。このように、入力信号Xの大きさの変化に応じて外乱信号の大きさを制御することで、入力信号が優ガウス性を有する場合には、元の外乱信号の性質がガウス性もしくは優ガウス性のいずれの場合であっても、f(x)で除算された外乱信号の性質は強い優ガウス性を有することになる。これにより、第二実施形態の未知伝達系推定技術は適用可能な範囲を拡張することができる。

［その他の変形例］
式（12）または式（14）において、分布の中心への集中度合（精度）を表すα, β等の係数は、未知伝達系のインパルス応答の振幅値や外乱信号の大きさを想定して、あらかじめ固定値として与えてもよい。また、フィルタ係数の推定値や観測された出力信号の大きさに基づいて時間とともに逐次修正してもよい。あるいは、ＥＭアルゴリズム等に基づいて反復学習により得られた値を逐次与えてもよい。

入力信号と出力信号を複数の周波数成分の信号に分割し、分割された周波数成分ごとの信号に対して、上述のいずれかの実施形態を適用し、分割された周波数成分ごとの伝達特性を推定してもよい。

この発明の未知伝達系推定技術によれば、優ガウス性の強い伝達特性（インパルス応答）を有する未知伝達系をより高速に推定することができる。実際、インパルス応答は、インパルス信号の応答波形であるから、インパルスが入力された直後に大きな振幅を取りながらも徐々に減衰する特性を有する場合が多く、このような応答は優ガウス性であるとみなせる場合が多い。また、外乱信号についても優ガウス性を仮定して導かれているため、音響エコーキャンセラにおけるダブルトークのように、外乱信号が優ガウス性を有する音声信号の場合などにおいても、精度高く未知の伝達特性を推定することが可能となる。その結果、エコーやハウリングの発生を抑えることができ、快適なハンズフリー通話を提供することが可能となる。

従来の未知伝達系推定技術は、未知伝達系のインパルス応答の振幅値および未知伝達系からの出力信号観測時に混入する外乱信号の振幅値がガウス分布に基づいて分布する場合においては適切であると解釈できるが、実際のインパルス応答や外乱信号は必ずしもガウス分布に基づいて生成されているとは限らない。この発明の未知伝達系推定技術は、特に、ガウス分布よりも零近傍の値を多く取ったり突発的に大きな値を取ったりするような優ガウス性の性質を、推定すべきインパルス応答および推定を阻害する外乱信号の双方に仮定し、出力信号が観測された条件のもとでのインパルス応答の振幅値が取り得る事後確率を定め、従来技術とは異なる新たな推定更新式を導くことで、実環境においても精度の高い推定を可能とした。

この発明は上述の実施形態に限定されるものではなく、この発明の趣旨を逸脱しない範囲で適宜変更が可能であることはいうまでもない。上記実施形態において説明した各種の処理は、記載の順に従って時系列に実行されるのみならず、処理を実行する装置の処理能力あるいは必要に応じて並列的にあるいは個別に実行されてもよい。

［プログラム、記録媒体］
上記実施形態で説明した各装置における各種の処理機能をコンピュータによって実現する場合、各装置が有すべき機能の処理内容はプログラムによって記述される。そして、このプログラムをコンピュータで実行することにより、上記各装置における各種の処理機能がコンピュータ上で実現される。

この処理内容を記述したプログラムは、コンピュータで読み取り可能な記録媒体に記録しておくことができる。コンピュータで読み取り可能な記録媒体としては、例えば、磁気記録装置、光ディスク、光磁気記録媒体、半導体メモリ等どのようなものでもよい。

また、このプログラムの流通は、例えば、そのプログラムを記録したDVD、CD-ROM等の可搬型記録媒体を販売、譲渡、貸与等することによって行う。さらに、このプログラムをサーバコンピュータの記憶装置に格納しておき、ネットワークを介して、サーバコンピュータから他のコンピュータにそのプログラムを転送することにより、このプログラムを流通させる構成としてもよい。

このようなプログラムを実行するコンピュータは、例えば、まず、可搬型記録媒体に記録されたプログラムもしくはサーバコンピュータから転送されたプログラムを、一旦、自己の記憶装置に格納する。そして、処理の実行時、このコンピュータは、自己の記録媒体に格納されたプログラムを読み取り、読み取ったプログラムに従った処理を実行する。また、このプログラムの別の実行形態として、コンピュータが可搬型記録媒体から直接プログラムを読み取り、そのプログラムに従った処理を実行することとしてもよく、さらに、このコンピュータにサーバコンピュータからプログラムが転送されるたびに、逐次、受け取ったプログラムに従った処理を実行することとしてもよい。また、サーバコンピュータから、このコンピュータへのプログラムの転送は行わず、その実行指示と結果取得のみによって処理機能を実現する、いわゆるASP（Application Service Provider）型のサービスによって、上述の処理を実行する構成としてもよい。なお、本形態におけるプログラムには、電子計算機による処理の用に供する情報であってプログラムに準ずるもの（コンピュータに対する直接の指令ではないがコンピュータの処理を規定する性質を有するデータ等）を含むものとする。

また、この形態では、コンピュータ上で所定のプログラムを実行させることにより、本装置を構成することとしたが、これらの処理内容の少なくとも一部をハードウェア的に実現することとしてもよい。

１００ 未知伝達系
２００、３００ 未知伝達系推定装置
２０１ フィルタ係数記憶部
２０２ 出力信号推定部
２０３ 残差信号算出部
２０４ 入力信号二乗値算出部
２０５、３０３ フィルタ期待値算出部
３０１ 入力信号絶対値算出部
３０２ 残差信号絶対値算出部

Claims (5)

1. 未知伝達系の伝達特性を表すフィルタ係数を記憶するフィルタ係数記憶部と、
   上記未知伝達系への入力信号の絶対値の系列を算出する入力信号絶対値算出部と、
   上記入力信号の系列と上記フィルタ係数との畳み込み演算により推定出力信号を算出する出力信号推定部と、
   上記入力信号を入力とする上記未知伝達系の出力を観測することで得られる観測出力信号と上記推定出力信号との残差信号を算出する残差信号算出部と、
   上記残差信号の絶対値を算出する残差信号絶対値算出部と、
   上記入力信号の系列と上記入力信号の絶対値の系列と上記残差信号と上記残差信号の絶対値とを用いて、上記フィルタ係数および上記観測出力信号に含まれうる外乱信号の振幅値が優ガウス性分布に従うものとして上記観測出力信号が観測された条件のもとでフィルタ係数が取り得る期待値を算出し、上記フィルタ係数を更新するフィルタ期待値算出部と、
   を含む未知伝達系推定装置。
2. 請求項１に記載の未知伝達系推定装置であって、
   ||・||₁をL1ノルムとし、・^Tを行列の転置とし、_Wを上記フィルタ係数の更新前の値を要素として持つベクトルとし、Wを上記フィルタ係数が更新後に取り得る値を要素として持つベクトルとし、yを上記観測出力信号とし、Xを上記入力信号の系列を要素として持つベクトルとし、αを上記外乱信号の精度とし、Bを上記フィルタ係数の精度を対角成分にもつ対角行列とし、
   上記フィルタ期待値算出部は、上記フィルタ係数が次式で定義される事後確率に従うものとした場合の期待値を算出するものである、
   

   

   
   未知伝達系推定装置。
3. 請求項１に記載の未知伝達系推定装置であって、
   ||・||₁をL1ノルムとし、・^Tを行列の転置とし、_Wを上記フィルタ係数の更新前の値を要素として持つベクトルとし、Wを上記フィルタ係数が更新後に取り得る値を要素として持つベクトルとし、yを上記観測出力信号とし、Xを上記入力信号の系列を要素として持つベクトルとし、αを上記外乱信号の精度とし、Bを上記フィルタ係数の精度を対角成分にもつ対角行列とし、f(X)を上記入力信号の大きさに依存して非負の値を取る関数とし、
   上記フィルタ期待値算出部は、上記フィルタ係数が次式で定義される事後確率に従うものとした場合の期待値を算出するものである、
   

   

   
   未知伝達系推定装置。
4. フィルタ係数記憶部が、未知伝達系の伝達特性を表すフィルタ係数を記憶しており、
   入力信号絶対値算出部が、上記未知伝達系への入力信号の絶対値の系列を算出する入力信号絶対値算出ステップと、
   出力信号推定部が、上記入力信号の系列と上記フィルタ係数との畳み込み演算により推定出力信号を算出する出力信号推定ステップと、
   残差信号算出部が、上記入力信号を入力とする上記未知伝達系の出力を観測することで得られる観測出力信号と上記推定出力信号との残差信号を算出する残差信号算出ステップと、
   残差信号絶対値算出部が、上記残差信号の絶対値を算出する残差信号絶対値算出ステップと、
   フィルタ期待値算出部が、上記入力信号の系列と上記入力信号の絶対値の系列と上記残差信号と上記残差信号の絶対値とを用いて、上記フィルタ係数および上記観測出力信号に含まれうる外乱信号の振幅値が優ガウス性分布に従うものとして上記観測出力信号が観測された条件のもとでフィルタ係数が取り得る期待値を算出し、上記フィルタ係数を更新するフィルタ期待値算出ステップと、
   を含む未知伝達系推定方法。
5. 請求項１から３のいずれかに記載の未知伝達系推定装置としてコンピュータを機能させるためのプログラム。

Images