JP2013179434A - 縦続接続型伝達系パラメータ推定方法、縦続接続型伝達系パラメータ推定装置、プログラム - Google Patents

Abstract

【課題】非線形系と線形系との縦続接続型伝達系のパラメータを推定する技術を提供する。
【解決手段】非線形系推定パラメータの更新値と線形系推定パラメータの更新値とが満足すべき近似誤差方程式における無数に存在する解のうち、モデル次数の違いを反映させた上で、非線形系推定パラメータの更新値の大きさ（二乗値）と線形系推定パラメータの更新値の大きさ（二乗値）の和が最も小さくなる値の組合せを算出することで、非線形系推定パラメータの更新値と線形系推定パラメータの更新値との間のバランスの取れた一意の解を得る。
【選択図】図２

Description

本発明は、電気信号、音、振動などが伝播する伝達系が、非線形な特性を有する系と線形な特性を有する系との縦続接続とみなせる場合に、当該伝達系の未知のパラメータを推定する技術に関する。

電気信号、音、振動などが伝播する未知の伝達系の特性を推定することは、当該伝達系への任意の入力に対する当該伝達系の出力を予測する上で有益である。伝達系を推定する問題は、有限のパラメータで表現可能なモデルを予め与え、当該パラメータの値を推定する問題として扱うことが多い。予め与えるモデルは、推定対象の伝達系の物理的素性（入力と出力の関係）に基づき選定される。例えば、入力と出力の関係が線形と見なせる場合、線形のFIR(Finite Impulse Response)フィルタやIIR(Infinite Impulse Response)フィルタが、離散時間領域における伝達系のモデルとして選ばれる。

未知の伝達系はこのような線形伝達系であるとは限らず、未知の伝達系における入力と出力の関係が非線形な特性を有する場合がある。例えば、未知の伝達系が線形特性と非線形特性の複合形態で表現されるような場合が考えられる。具体的には、図１に示すように、推定対象の伝達系１が非線形系２と線形系３の縦続接続で表現できる場合であって、伝達系１への入力は非線形系２に入力され、非線形系２の出力は線形系３に入力され、線形系３の出力が伝達系１の出力となっている。

このような伝達系１において、それぞれ独立した系の非線形系２と線形系３とによる縦続接続型のモデルとして与えることができ、なおかつ、非線形系２の特性は過去の入力（出力）に依存せず、線形系３の特性は有限長のインパルス応答と見なせる場合における、伝達系１のパラメータ推定の問題を、非特許文献１に開示される技術は扱っている。

[非線形系のモデル]
図１の非線形系２への入力信号をx(n)、非線形系２の出力信号をs(n)とすると、式（１）なる関係が与えられる。ここで、nは離散時間、fは非線形関数であり、出力s(n)は、入力x(n)に対しては、同じ離散時間nにおける値にのみ依存するとして関係づけられている。
s(n)=f(x(n)) （１）

ここで、非線形関数fの特性が未知である場合に、その特性を直接的に推定することは困難であるため、従来から、非線形関数fを有限数のモデルパラメータで近似した上で、パラメータ推定する手法がある。例えば、P次の級数展開の有限打ち切り近似により、式（２）で近似した上で、モデルパラメータa₀，a₁，a₂，…，a_Pを推定するという手法である。ここで、Pは0以上の整数で、関数g₀，g₁，g₂，…，g_Pは、互いに異なり、なおかつ、非線形関数fの特性とは独立して与えられているとする。

非特許文献１には、関数g₀，g₁，g₂，…，g_Pの例としてべき級数展開の有限項の打ち切りに相当する式（３）が示されている。
g_p(x(n))=x(n)^p, p=0,1,2,…,P （３）

[線形系のモデル]
図１の線形系３への入力信号をs(n)とその出力信号y(n)とすると、y(n)は式（４）の畳込演算の形で表すことができる。ここで有限長のインパルス応答h₀，h₁，h₂，…，h_Lが推定すべき未知のモデルパラメータである。ただし、次数Lは0以上の整数であり、s(n)は、実際には直接観測できない内部信号である。

以上により、推定対象のパラメータは、式（２）におけるa₀，a₁，a₂，…，a_Pと式（４）におけるh₀，h₁，h₂，…，h_Lである。

非特許文献１には、このような枠組におけるパラメータ推定方法の一例が示されている。非特許文献１の方法によると、離散時間n+1で用いる級数展開係数（非線形系パラメータ）a₀，a₁，a₂，…，a_Pの推定値（非線形系推定パラメータ）a^₀(n)，a^₁(n)，a^₂(n)，…，a^_P(n)は式（５）に基づいて逐次更新によって推定される。

同様にインパルス応答（線形系パラメータ）h₀，h₁，h₂，…，h_Lの推定値（線形系推定パラメータ）h^₀(n)，h^₁(n)，h^₂(n)，…，h^_L(n)は式（６）に基づいて逐次更新によって推定される。

ここで、A(n)，H(n)，X(n)，S(n)，G(n)はそれぞれ式（７）−（１１）で表される。Tはベクトルまたは行列の転置、μ_a，μ_hは、それぞれパラメータの更新量を調整するステップサイズである。‖・‖はノルムを表す。

また、e(n)は、式（１２）で表される誤差信号である。

A. Stenger, W. Kellermann; "Nonlinear acoustic echo cancellation with fast converging memoryless preprocessor", IEEE International Conference on Acoustics, Speech, and Signal Processing, 2000. (ICASSP '00), Vol. 2, II805 - II808, 2000.

図１に示すような非線形系と線形系とが縦続接続された伝達系のパラメータを、式（５）および式（６）に従って逐次更新して推定する従来技術では、正規化（‖X(n)‖²または‖S(n)‖²による除算）が非線形系推定パラメータa^₀(n)，a^₁(n)，a^₂(n)，…，a^_P(n)と、線形系推定パラメータの推定値h^₀(n)，h^₁(n)，h^₂(n)，…，h^_L(n)に対して、独立に行われており、相互のバランスが考慮されていない。このため、予め人手で選定して与える必要のある非線形系パラメータ用のステップサイズμ_aおよび線形系パラメータ用のステップサイズμ_hの間の関係が明確でなく、μ_aとμ_hの適切な値の組合せを簡単に選び出すことができなかった。また、非線形系推定パラメータの次数Pの大きさと線形系推定パラメータの次数Lの大きさとの違いも考慮されていない。

そこで本発明は、非線形系推定パラメータの次数Pの大きさと線形系推定パラメータの次数Lの大きさとの違いを考慮した推定パラメータの逐次更新に基づいて、非線形系と線形系との縦続接続型伝達系のパラメータを推定する技術を提供することを目的とする。

本発明は、過去の入出力に依存しない特性を持つ非線形系と有限長のインパルス応答特性を持つ線形系との縦続接続とみなせる伝達系について、当該非線形系のモデルパラメータおよび当該線形系のモデルパラメータを推定する縦続接続型伝達系パラメータ推定技術であって、Pをゼロ以上の整数として、時刻nにおける伝達系への入力信号x(n)に対し、P+1個の異なる関数g₀，g₁，g₂，…，g_Pを適用した結果と、当該各関数に対応する非線形系のモデルパラメータの推定値a^₀(n)，a^₁(n)，a^₂(n)，…，a^_P(n)との積和として、非線形系の出力推定値s^(n)を得る処理と、Lをゼロ以上の整数として、時刻n-Lから時刻nまでの各時刻に対応する非線形系の出力推定値s^(n-L)，…，s^(n)と、当該各出力推定値に対応する線形系のモデルパラメータの推定値h^₀(n)，h^₁(n)，h^₂(n)，…，h^_L(n)との積和として、線形系の出力推定値y^(n)を得る処理と、伝達系の出力信号y(n)から線形系の出力推定値y^(n)を引いた誤差信号e(n)を得る処理と、時刻n-Lから時刻nまでの各時刻に対応する各入力信号x(n-L)，…，x(n)を、P+1個の各関数g₀，g₁，g₂，…，g_Pに適用し、当該各関数毎に得られるL+1個の結果と、線形系のL+1個のモデルパラメータの推定値h^₀(n)，h^₁(n)，h^₂(n)，…，h^_L(n)との積和に相当する値を時間積和値として得、当該各関数毎のP+1個の時間積和値を要素として持つ非線形系入力評価ベクトルX(n)と、時刻n-Lから時刻nまでの各時刻に対応する各入力信号x(n-L)，…，x(n)を、P+1個の各関数g₀，g₁，g₂，…，g_Pに適用し、時刻n-Lから時刻n毎に得られるP+1個の結果と、非線形系のモデルパラメータの推定値a^₀(n)，a^₁(n)，a^₂(n)，…，a^_P(n)との積和に相当する値を関数積和値として得、当該時刻毎のL+1個の関数積和値を要素として持つ線形系入力評価ベクトルS(n)と、誤差信号e(n)、非線形系のモデル次数P、線形系のモデル次数Lを用い、Tをベクトルの転置として、非線形系のモデル次数Pおよび線形系のモデル次数Lに応じて、PがLよりも小さいほど更新量が大きくなるように正規化された非線形系のモデルパラメータの更新量で、非線形系入力評価ベクトルX(n)をスカラー倍したベクトルを非線形系のモデルパラメータベクトル[a^₀(n)，a^₁(n)，a^₂(n)，…，a^_P(n)]^Tに加えることで、時刻n+1における非線形系のモデルパラメータベクトル[a^₀(n+1)，a^₁(n+1)，a^₂(n+1)，…，a^_P(n+1)]^Tを推定し、非線形系のモデル次数Pおよび線形系のモデル次数Lに応じて、PがLよりも小さいほど更新量が小さくなるように正規化された線形系のモデルパラメータの更新量で、線形系入力評価ベクトルS(n)をスカラー倍したベクトルを線形系のモデルパラメータベクトル[h^₀(n)，h^₁(n)，h^₂(n)，…，h^_L(n)]^Tに加えることで、時刻n+1における線形系のモデルパラメータベクトル[h^₀(n+1)，h^₁(n+1)，h^₂(n+1)，…，h^_L(n+1)]^Tを推定する処理を行う。

本発明に拠ると、詳しくは後の説明に譲るが、非線形系推定パラメータの更新量と線形系推定パラメータの更新量との相互のバランスをモデル次数に応じて正規化しており、予め人手で選定して与える必要のある更新量調整のためのステップサイズを非線形系と線形系について個別に選定することなく等しい値を与えても推定値の収束速度が非線形系と線形系とで偏りがなく、適切な性能を達成することが可能となる。

伝達系の構成を示す図。 縦続接続型伝達系パラメータ推定装置１００の機能構成例を示す図。 縦続接続型伝達系パラメータ推定装置１００をエコーキャンセラとして機能させる場合の機能構成例を示す図。

[本発明の対象となる未知伝達系]
本発明の対象となる未知伝達系（推定対象の伝達系とも云う）は、図１に示す伝達系１のように、「それぞれ独立した系である非線形系２と線形系３の縦続接続で伝達系１が表現されており、伝達系１への入力は非線形系２に入力され、非線形系２の出力は線形系３に入力され、線形系３の出力が伝達系１の出力となっており、かつ、非線形系２の特性は過去の入力（出力）に依存せず、線形系３の特性は有限長のインパルス応答と見なせる」伝達系である。

[非線形系のモデル]
既述のように、図１の非線形系２への入力信号をx(n)、非線形系２の出力信号をs(n)とすると、上記式（１）なる関係が与えられる。記号も踏襲し、nは離散時間、fは非線形関数を表すとする。出力s(n)は、入力x(n)に対しては、同じ離散時間nにおける値にのみ依存するとして関係づけられている。

また、既述のように、非線形関数fを有限数のモデルパラメータで近似する。ここでは、P次の級数展開の有限打ち切り近似により、上記式（２）で近似する。記号も踏襲し、Pは0以上の整数で、関数g₀，g₁，g₂，…，g_Pは、互いに異なり、なおかつ、非線形関数fの特性とは独立して与えられているとする。

[線形系のモデル]
既述のように、図１の線形系３への入力信号をs(n)とその出力信号y(n)とすると、y(n)は上記式（４）の畳込演算の形で表すことができる。記号も踏襲し、次数Lは0以上の整数であり、s(n)は、実際には直接観測できない内部信号として扱われる。

以上により、伝達系１にて推定対象となるモデルパラメータは、上記式（２）における関数g₀，g₁，g₂，…，g_Pの係数a₀，a₁，a₂，…，a_Pと上記式（４）における有限長のインパルス応答h₀，h₁，h₂，…，h_Lである。

また、ここでの説明ではA(n)，H(n)，X(n)，S(n)，G(n)はそれぞれ上記式（７）−（１１）で表されるとする。記号も踏襲し、Tはベクトルまたは行列の転置である。‖・‖はノルムを表す。また、e(n)を上記式（１２）で表される誤差信号とする。

さて、非線形系推定パラメータa^₀(n)，a^₁(n)，a^₂(n)，…，a^_P(n)と線形系推定パラメータh^₀(n)，h^₁(n)，h^₂(n)，…，h^_L(n)とを更新する際、両者の相対的な更新量の大きさを規定することを考える。

前提として、式（１２）の誤差信号について式（１３）で与えられる評価式について、バラメータベクトルA(n)、H(n)による偏微分式（１４），（１５）の符号を正負入れ替えたものが、パラメータベクトルの更新方向を与える。

しかしながら、これらの方向にどれだけの大きさで更新すればよいのかは、この時点では明確ではない。そこで、パラメータベクトルの更新量の適切な大きさを議論するために、従来技術における式（５），（６）に相当するパラメータ更新式を式（１６），（１７）で与える。

ここでΔA(n)，ΔH(n)は、各パラメータベクトルの更新ベクトルで、これらの方向は、それぞれ、式（１４），（１５）のベクトルと等しい（符号が反対の場合も含む）とすれば、α、βをスカラー値として、式（１８），（１９）のように表すことができる。

ここで、更新した後のパラメータベクトルA(n+1)、H(n+1)によって、時刻nの出力信号が正しく模擬されるためには式（２０）を満足する必要がある。
y(n)=H(n+1)^TG(n)A(n+1) （２０）

式（２０）に式（１２），（１６），（１７），（１８），（１９）を適用すると、誤差方程式（２１）が得られる。

ここで、更新の大きさは小さいと仮定して、式（２１）の右辺第三項を無視すると、近似誤差方程式として式（２２）が得られる。

式（２２）を満足するα、βの組合せは無数に存在するが、その中で、α²＋β²を最小にするα、βは一意に存在して、αは式（２４）で、βは式（２５）で与えられる。

これらを式（１８），（１９）へ、それぞれ代入すると、更新ベクトルは式（２６），（２７）で表される。

[パラメータの次数に基づく正規化の厳密化]
このままでは、非線形系パラメータPと線形系パラメータの次数Lとの違いが考慮されていないため、次数の低いパラメータのグループの更新量が小さくなり、推定すべき特性への収束が遅くなるという問題が残る。多くの場合、非線形系のモデルは少ない次数で表現されることが多いが、線形系のモデルは、室内の音響インパルス応答などに相当し、次数が非常に高くなることが多い。そこで、Lの値が数百から数万、Pの値1から数百程度でLがPよりも10〜数100倍程度大きい場合でも推定性能が劣化しないように考慮する。そこで、本発明では、式（１８），（１９）の代わりに、更新ベクトルを式（２８），（２９）で与える。

このとき、式（２４），（２５）を得る議論と同様の議論を経ると、αは式（３０）で、βは式（３１）で表される。

このとき、更新ベクトルは、式（３２），（３３）で表される。

これらの更新ベクトルは、非線形系パラメータと線形系パラメータの次数の違いも考慮に入れて正規化されているため、実際の更新計算において、個別に異なるステップサイズを与えずにパラメータの推定が可能となる。すなわち、共通のステップサイズμを導入し、式（３４），（３５）として、更新することが可能となる。

[パラメータ次数に基づく正規化の一般的な実施形態]
式（３４），（３５）の本質は、非線形系、線形系それぞれのパラメータ次数の違いがあっても、パラメータ更新に偏りが起きないように正規化することであり、式（３４），（３５）と等価、あるいは同等の更新式も、本発明に包含される。

例えば、式（３６），（３７）を例示できる。w(P,L)は、PがLよりも大きい場合にその大きさの違いに応じて１より大きな正の値をとり、PとLがほぼ等しい場合に１付近の正の値をとり、PがLよりも小さい場合にその小ささの違いに応じて1より小さい正の値をとるような任意の関数である。w(P,L)の値の設定を、予め用意されたルックアップテーブルを参照することによって、あるいは、条件分岐などによって実現してもよい。w(L,P)は、w(P,L)におけるPとLの関係を入れ替えて同様に実現する。

また、一般式としての式（３６），（３７）における‖X(n)‖²＋w(P,L)‖S(n)‖²または、‖S(n)‖²＋w(L,P)‖X(n)‖²による除算は、それらの値が0となった場合に除算結果が発散するのを防ぐため、小さな正の定数δを導入し、‖X(n)‖²＋w(P,L)‖S(n)‖²＋δや‖S(n)‖²＋w(L,P)‖X(n)‖²＋δによる除算で置き換えても、本発明の効果は維持される。

[推定パラメータの一部を定数とみなす実施形態]
さらに、本発明では、一部の推定パラメータを定数とみなし固定することで、推定パラメータの一意な解への収束性を高めることができる。

例えば、非線形系推定パラメータのうち、a^₀(n)=1に固定すると、式（３６）などのA(n)の更新式において、A(n)を式（３８）に置き換え、更新式に現れるX(n)についても式（３９）に置き換えた上で、式（９）により計算することになる。

このとき、A(n)の次数が一つ減少したと見なせるので、式（３４）の更新式を式（４０）に置き換えるのが、より厳密である。

他方、このときのH(n)の更新についても、式（３５）の更新式を式（４１）に置き換えるのが、より厳密である。

ただし、式（４０），（４１）で用いるS(n)は、式（３８)，（３９）のA(n)、G(n)を用いて式（４２）で計算することに相当する。

[非線形系の級数展開関数の例]
本発明では、非線形系と線形系のパラメータの次数の違いに着目することで伝達系の推定精度を高めている。すなわち、式（２）の近似式における関数g₀，g₁，g₂，…，g_Pの与え方に依らず、従来の更新式と比較して優位となる効果は変わらない。

具体的な関数g₀，g₁，g₂，…，g_Pとして、上記式（３）のべき級数展開を用いてもよいし、そのほか、フーリエ級数や式（４３）のような振幅値の大きさで信号を分割する関数を用いてもよい。ただし、p=0,1,2,....,P、q=0,1,2,...,Pである。振幅分割閾値は、0=T₀<T₁<T₂<…<T_P+1であり、例えば、16ビットのPCMデジタル信号を対象とする場合は、上限値T_P+1を32768に相当させれば全ての信号振幅値を網羅できる。

p=qとして、関数gの添え字pと閾値Tの添え字qの順序を揃えてもよいが、一般的には、1からPの値が一度だけ選ばれるようにしさえすれば、pとqの値の与え方は任意でよい。また、式（４３）と等価な出力を与える他の関数表現や条件分岐、ルックアップテーブル参照等によって関数g₀，g₁，g₂，…，g_Pを実現してもよい。

また、振幅値の符号に対して非対称な分割として、式（４４），（４５）のような関数を用いてもよい。ただし、k=0,1,2,...,K、i=0,1,2,...,K、m=0,1,2,...,M、r=0,1,2,...,Mであり、K+M=P-1の関係のもとで、次数Pが与えられる。なお、KとMの値は必ずしも等しい必要はなく、閾値は0=V₀<V₁<V₂<…<V_K+1，0=U₀<U₁<U₂<…<U_M+1の関係を満足すればよい。

さらに、振幅値の符号による対称な分割と非対称な分割を組合わせてもよく、信号x(n)の絶対値が、ある大きさの範囲において、式（４３）に基づき符号に対して対称な分割を行い、それ以外の範囲においては、式（４４），（４５）のように符号に対して非対称な分割を行ってもよい。

[縦続接続型伝達系パラメータ推定装置／方法]
以上の説明を踏まえて、一例として図３に示すような縦続接続型伝達系パラメータ推定装置１００を構成することができる。

現在の離散時間nにおいて、推定対象の伝達系１への入力信号x(n)が非線形関数級数展開部５０１に入力され、非線形関数級数展開部５０１は、P次の展開関数g₀，g₁，g₂，…，g_Pを用いて、P+1個の級数展開信号g₀(x(n))，g₁(x(n))，g₂(x(n))，…，g_P(x(n))を得る。ここで、展開関数は、上記のべき級数展開関数、フーリエ級数展開関数、閾値分割関数などが適用される（式（３），（４３），（４４），（４５）参照）。

非線形パラメータ乗算・合成部５０２は、離散時間nに関するP+1個の信号g₀(x(n))，g₁(x(n))，g₂(x(n))，…，g_P(x(n))と離散時刻nにおいて推定されているP+1個の非線形系推定パラメータa^₀(n)，a^₁(n)，a^₂(n)，…，a^_P(n)との積和、つまり推定対象の伝達系１に含まれる非線形系２の出力に相当する推定値s^(n)=a^₀(n)g₀(x(n))+a^₁(n)g₁(x(n))+a^₂(n)g₂(x(n))+…+a^_P(n)g_P(x(n))を出力する。

この結果、信号蓄積部５０３には、離散時間n-Lから離散時間nまでのL+1個の非線形系出力推定値s^(n),s^(n-1),s^(n-2),…,s^(n-L)が記憶されている。

線形パラメータ乗算・合成部５０４は、信号蓄積部５０３から取得したL+1個の非線形系出力推定値s^(n),s^(n-1),s^(n-2),…,s^(n-L)と離散時刻nにおいて推定されているL+1個の線形系推定パラメータh^₀(n)，h^₁(n)，h^₂(n)，…，h^_L(n)との積和、つまり推定対象の伝達系１からの出力信号の推定値y^(n)=h^₀(n)s^(n)+h^₁(n)s^(n-1)+h^₂(n)s^(n-2)+…+h^_L(n)s^(n-L)を出力する。

誤差算出部５０５は、推定対象の伝達系１からの出力信号y(n)から、その推定値y^(n)を差し引き、誤差信号e(n)を出力する。

さらに、非線形・線形パラメータ統合正規化更新部５０７は、非線形関数級数展開部５０１から離散時刻nごとに級数展開信号g₀(x(n))，g₁(x(n))，g₂(x(n))，…，g_P(x(n))を取得しており、離散時刻n-Lから離散時刻nまでの、(P+1)×(L+1)個の蓄積級数展開信号g₀(x(n))，g₁(x(n))，g₂(x(n))，…，g_P(x(n))，g₀(x(n-1))，g₁(x(n-1))，g₂(x(n-1))，…，g_P(x(n-1))，…，g₀(x(n-L))，g₁(x(n-L))，g₂(x(n-L))，…，g_P(x(n-L))を蓄積している。そして、非線形・線形パラメータ統合正規化更新部５０７は、離散時刻nにおいて推定されているL+1個の線形系推定パラメータh^₀(n)，h^₁(n)，h^₂(n)，…，h^_L(n)と各次数p(p=0,1,2,....,P)に対応したL+1個の蓄積級数展開信号g_p(x(n))，g_p(x(n-1))，g_p(x(n-2))，…，g_p(x(n-L))との積和演算結果を要素として持つ(P+1)次元の非線形系入力評価ベクトルX(n)を得る。

さらに、非線形・線形パラメータ統合正規化更新部５０７は、離散時刻nにおいて推定されているP+1個の非線形系推定パラメータa^₀(n)，a^₁(n)，a^₂(n)，…，a^_P(n)と、離散時刻n-j(j=0,1,2,...,L)に対応したP+1個の蓄積級数展開信号g₀(x(n-j))，g₁(x(n-j))，g₂(x(n-j))，…，g_P(x(n-j))との積和演算結果を要素として持つ(L+1)次元の線形系入力評価ベクトルS(n)を直接計算する。

あるいは、非線形・線形パラメータ統合正規化更新部５０７は、信号蓄積部５０３に記憶されている、離散時間n-Lから離散時間nまでのL+1個の非線形系出力推定値s^(n),s^(n-1),s^(n-2),…,s^(n-L)を用いて、線形系入力評価ベクトルS(n)を得ることもできる。

そして、非線形・線形パラメータ統合正規化更新部５０７は、非線形系入力評価ベクトルX(n)、線形系入力評価ベクトルS(n)、誤差信号e(n)、非線形系モデルの次数P、線形系モデルの次数Lを用いて、非線形系推定パラメータを式（４６）により更新し、線形系推定パラメータを式（４７）により更新する。

これらの更新結果は伝達系１の特性の推定結果であり、縦続接続型伝達系パラメータ推定装置１００は必要に応じて当該推定結果を出力する。更新結果の出力は、離散時間ｎの更新毎に遂次実行してもよいし、伝達系１の時間変動が生じない範囲で一定更新回数（例えばＬ×数１０回）経過後に実行してもよい。また、縦続接続型伝達系パラメータ推定装置１００は、非線形系パラメータ、あるいは線形系パラメータのいずれか一方のみを出力してもよい。

なお、式（４６）においてw(P,L)は、既述のように、例えば、PがLよりも大きい場合にその大きさの違いに応じて１より大きな正の値をとり、PとLがほぼ等しい場合に１付近の正の値をとり、PがLよりも小さい場合にその小ささの違いに応じて1より小さい正の値をとるような任意の関数である。また、式（４７）においてw(L,P)は、既述のように、w(P,L)におけるPとLの関係を入れ替えて同様に実現する。また、δは、零除算防止のための正の定数であるが、必ずしも必要となるものではない。

さらに、μ_a，μ_hは更新量を調整する0から１の範囲の定数である。μ_a＝μ_hとしてもよいし、μ_a＝μ_hのときの収束特性を基準として、このときのμ_aまたはμ_hの値を微調整したμ'_aまたはμ'_hを用いてもよい。

式（４６）で得られた非線形系推定パラメータの更新値a^₀(n+1)，a^₁(n+1)，a^₂(n+1)，…，a^_P(n+1)は、離散時間n+1において、非線形パラメータ乗算・合成部５０２で用いられる。式（４７）で得られた線形系推定パラメータの更新値h^₀(n+1)，h^₁(n+1)，h^₂(n+1)，…，h^_L(n+1)は、離散時間n+1において、線形パラメータ乗算・合成部５０４で用いられる。

もちろん、式（４６），（４７）に替えて、式（３４），（３５）や式（４０），（４１）を用いてもよい。

このように、本発明の縦続接続型伝達系パラメータ推定技術によると、非線形系推定パラメータの更新値と線形系推定パラメータの更新値とが満足すべき近似誤差方程式における無数に存在する解のうち、モデル次数の違いを反映させた上で、非線形系推定パラメータの更新値の大きさ（二乗値）と線形系推定パラメータの更新値の大きさ（二乗値）の和が最も小さくなる値の組合せを算出することで、非線形系推定パラメータの更新値と線形系推定パラメータの更新値との間のバランスの取れた一意の解を得ることができる。

［音響エコーキャンセラへの適用］
本発明の縦続接続型伝達系パラメータ推定技術は、例えば音響エコーキャンセラへの転用が可能である。音響エコーキャンセラは、スピーカから再生される音響信号が室内で反響し、同一室内にあるマイクロホンへ回り込むエコーを推定して得られる擬似エコー信号をマイクロホンの収音信号から差し引くことで、エコーの消去された信号を出力するものである。図３に示すように、スピーカで再生するための再生信号x(n)に増幅器２０１でゲインをかけてスピーカ２０２から再生する際、増幅器２０１やスピーカ２０２により信号に歪が生じる場合がある。すなわち、スピーカで再生するための再生信号x(n)が、歪を発生させる増幅器２０１やスピーカ２０２を経由し、スピーカから放音される信号が得られるまでの伝達特性は非線形とみなし、非線形系２を構成しているものと想定できる。また、スピーカ２０２から放出された音響信号が、室内の反響路３０１を経由して、マイクロホン３０２に回り込み、収音信号y(n)が得られるまでの伝達特性は線形とみなし、線形系３を構成するものと想定できる。

この場合、縦続接続型伝達系パラメータ推定装置１００は音響エコーキャンセラ１００ａとして動作することが可能であり、上述の説明にて入力信号x(n)をスピーカで再生するための再生信号x(n)に、出力信号y(n)をスピーカと同一の空間で収音した収音信号y(n)に読み替えればよいので、重複説明を省略する。なお、音響エコーキャンセラ１００ａは、非線形系推定パラメータも線形系推定パラメータも外部に出力する必要はない代わりに、誤差信号e(n)を送信信号として出力する。

収音信号y(n)には、伝達系に由来するエコーのほかに、近端の話者の発話音声などが含まれる。近端の話者の音声とエコーが混合している場合には、y(n)から擬似エコー信号である推定対象の伝達系１からの出力信号の推定値y^(n)を差し引いて得られる誤差信号e(n)は、近端の話者の音声が主となり、通信相手への送信信号となる。しかしながら同時に、誤差信号e(n)の中に存在する近端話者の音声は、非線形系パラメータや線形系パラメータの推定結果を乱す要因ともなる。そこで、ステップサイズμの値を固定値とする場合は、その値を十分に小さく与えたり、近端話者の音声の存在を検知できる場合は、近端話者の音声が検知された場合にのみ、ステップサイズμを小さくするように制御したりする対処が別途必要である。

＜縦続接続型伝達系パラメータ推定装置のハードウェア構成例＞
上述の実施形態に関わる縦続接続型伝達系パラメータ推定装置は、キーボードなどが接続可能な入力部、液晶ディスプレイなどが接続可能な出力部、ＣＰＵ（Central Processing Unit）〔キャッシュメモリなどを備えていてもよい。〕、メモリであるＲＡＭ（Random Access Memory）やＲＯＭ（Read Only Memory）と、ハードディスクである外部記憶装置、並びにこれらの入力部、出力部、ＣＰＵ、ＲＡＭ、ＲＯＭ、外部記憶装置間のデータのやり取りが可能なように接続するバスなどを備えている。また必要に応じて、縦続接続型伝達系パラメータ推定装置に、ＣＤ−ＲＯＭなどの記憶媒体を読み書きできる装置（ドライブ）などを設けるとしてもよい。このようなハードウェア資源を備えた物理的実体としては、汎用コンピュータなどがある。

縦続接続型伝達系パラメータ推定装置の外部記憶装置には、パラメータを推定するためのプログラム並びにこのプログラムの処理において必要となるデータなどが記憶されている〔外部記憶装置に限らず、例えばプログラムを読み出し専用記憶装置であるＲＯＭに記憶させておくなどでもよい。〕。また、これらのプログラムの処理によって得られるデータなどは、ＲＡＭや外部記憶装置などに適宜に記憶される。以下、データやその格納領域のアドレスなどを記憶する記憶装置を単に「記憶部」と呼ぶことにする。

縦続接続型伝達系パラメータ推定装置の記憶部には、入力信号x(n)を用いて級数展開信号を得るためのプログラム、級数展開信号と非線形系推定パラメータとの積和（推定対象の伝達系に含まれる非線形系の出力に相当する推定値）を求めるためのプログラム、非線形系出力推定値と線形系推定パラメータとの積和（推定対象の伝達系からの出力信号の推定値）を求めるためのプログラム、推定対象の伝達系からの出力信号から、その推定値を差し引き、誤差信号を出力するためのプログラム、非線形系入力評価ベクトル、線形系入力評価ベクトル、誤差信号、非線形系モデルの次数、線形系モデルの次数を用いて、非線形系推定パラメータと線形系推定パラメータを更新するためのプログラムが記憶されている。

縦続接続型伝達系パラメータ推定装置では、記憶部に記憶された各プログラムとこの各プログラムの処理に必要なデータが必要に応じてＲＡＭに読み込まれて、ＣＰＵで解釈実行・処理される。この結果、ＣＰＵが所定の機能（非線形関数級数展開部、非線形パラメータ乗算・合成部、線形パラメータ乗算・合成部、誤差算出部、非線形・線形パラメータ統合正規化更新部）を実現することでパラメータ推定が実現される。

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

また、上記実施形態において説明したハードウェアエンティティ（縦続接続型伝達系パラメータ推定装置）における処理機能をコンピュータによって実現する場合、ハードウェアエンティティが有すべき機能の処理内容はプログラムによって記述される。そして、このプログラムをコンピュータで実行することにより、上記ハードウェアエンティティにおける処理機能がコンピュータ上で実現される。

この処理内容を記述したプログラムは、コンピュータで読み取り可能な記録媒体に記録しておくことができる。コンピュータで読み取り可能な記録媒体としては、例えば、磁気記録装置、光ディスク、光磁気記録媒体、半導体メモリ等どのようなものでもよい。具体的には、例えば、磁気記録装置として、ハードディスク装置、フレキシブルディスク、磁気テープ等を、光ディスクとして、ＤＶＤ（Digital Versatile Disc）、ＤＶＤ−ＲＡＭ（Random Access Memory）、ＣＤ−ＲＯＭ（Compact Disc Read Only Memory）、ＣＤ−Ｒ（Recordable）／ＲＷ（ReWritable）等を、光磁気記録媒体として、ＭＯ（Magneto-Optical disc）等を、半導体メモリとしてＥＥＰ−ＲＯＭ（Electronically Erasable and Programmable-Read Only Memory）等を用いることができる。

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

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

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

Claims (9)

1. 過去の入出力に依存しない特性を持つ非線形系と有限長のインパルス応答特性を持つ線形系との縦続接続とみなせる伝達系について、当該非線形系のモデルパラメータおよび当該線形系のモデルパラメータを推定する縦続接続型伝達系パラメータ推定方法であって、
   Pをゼロ以上の整数として、時刻nにおける上記伝達系への入力信号x(n)に対し、P+1個の異なる関数g₀，g₁，g₂，…，g_Pを適用した結果と、当該各関数に対応する上記非線形系のモデルパラメータの推定値a^₀(n)，a^₁(n)，a^₂(n)，…，a^_P(n)との積和として、上記非線形系の出力推定値s^(n)を得るステップと、
   Lをゼロ以上の整数として、時刻n-Lから時刻nまでの各時刻に対応する上記非線形系の出力推定値s^(n-L)，…，s^(n)と、当該各出力推定値に対応する上記線形系のモデルパラメータの推定値h^₀(n)，h^₁(n)，h^₂(n)，…，h^_L(n)との積和として、上記線形系の出力推定値y^(n)を得るステップと、
   上記伝達系の出力信号y(n)から上記線形系の出力推定値y^(n)を引いた誤差信号e(n)を得るステップと、
   時刻n-Lから時刻nまでの各時刻に対応する各入力信号x(n-L)，…，x(n)を、P+1個の上記各関数g₀，g₁，g₂，…，g_Pに適用し、当該各関数毎に得られるL+1個の結果と、上記線形系のL+1個のモデルパラメータの推定値h^₀(n)，h^₁(n)，h^₂(n)，…，h^_L(n)との積和に相当する値を時間積和値として得、当該各関数毎のP+1個の上記時間積和値を要素として持つ非線形系入力評価ベクトルX(n)と、時刻n-Lから時刻nまでの各時刻に対応する各入力信号x(n-L)，…，x(n)を、P+1個の上記各関数g₀，g₁，g₂，…，g_Pに適用し、上記時刻n-Lから上記時刻n毎に得られるP+1個の結果と、上記非線形系のモデルパラメータの推定値a^₀(n)，a^₁(n)，a^₂(n)，…，a^_P(n)との積和に相当する値を関数積和値として得、当該時刻毎のL+1個の上記関数積和値を要素として持つ線形系入力評価ベクトルS(n)と、上記誤差信号e(n)、上記非線形系のモデル次数P、上記線形系のモデル次数Lを用い、Tをベクトルの転置として、上記非線形系のモデル次数Pおよび上記線形系のモデル次数Lに応じて、PがLよりも小さいほど更新量が大きくなるように正規化された上記非線形系のモデルパラメータの更新量で、上記非線形系入力評価ベクトルX(n)をスカラー倍したベクトルを上記非線形系のモデルパラメータベクトル[a^₀(n)，a^₁(n)，a^₂(n)，…，a^_P(n)]^Tに加えることで、時刻n+1における上記非線形系のモデルパラメータベクトル[a^₀(n+1)，a^₁(n+1)，a^₂(n+1)，…，a^_P(n+1)]^Tを推定し、上記非線形系のモデル次数Pおよび上記線形系のモデル次数Lに応じて、PがLよりも小さいほど更新量が小さくなるように正規化された上記線形系のモデルパラメータの更新量で、上記線形系入力評価ベクトルS(n)をスカラー倍したベクトルを上記線形系のモデルパラメータベクトル[h^₀(n)，h^₁(n)，h^₂(n)，…，h^_L(n)]^Tに加えることで、時刻n+1における上記線形系のモデルパラメータベクトル[h^₀(n+1)，h^₁(n+1)，h^₂(n+1)，…，h^_L(n+1)]^Tを推定する更新ステップと
   を有する縦続接続型伝達系パラメータ推定方法。
2. 請求項１に記載の縦続接続型伝達系パラメータ推定方法であって、
   上記更新ステップにおいて、μをパラメータの更新量を調整するステップサイズ、‖・‖はノルムを表すとして、
   

   

   
   であることを特徴とする縦続接続型伝達系パラメータ推定方法。
3. 請求項１に記載の縦続接続型伝達系パラメータ推定方法であって、
   上記更新ステップにおいて、μをパラメータの更新量を調整するステップサイズ、‖・‖はノルムを表すとして、
   

   

   
   、ただし、w(P,L)およびw(P,L)はPとLとの相違に基づいて正値をとる
   ことを特徴とする縦続接続型伝達系パラメータ推定方法。
4. 請求項１に記載の縦続接続型伝達系パラメータ推定方法であって、
   上記更新ステップにおいて、a^₀(n)=1として固定値で与え、上記X(n)の上記関数g₀に関する要素を除いたP次元のベクトルとし、μをパラメータの更新量を調整するステップサイズ、‖・‖はノルムを表すとして、
   

   

   
   であることを特徴とする縦続接続型伝達系パラメータ推定方法。
5. 請求項１から請求項４のいずれかに記載の縦続接続型伝達系パラメータ推定方法であって、
   上記入力信号はスピーカで再生するための再生信号であり、
   上記出力信号は上記スピーカと同一の空間で収音された収音信号であり、
   上記誤差信号は上記収音信号からエコーが除去された送話信号である
   ことを特徴とする縦続接続型伝達系パラメータ推定方法。
6. 請求項１から請求項５のいずれかに記載の縦続接続型伝達系パラメータ推定方法であって、
   上記関数g_pは、
   

   

   
   、ただし、p=0,1,2,....,P、q=0,1,2,...,Pである
   ことを特徴とする縦続接続型伝達系パラメータ推定方法。
7. 請求項１から請求項５のいずれかに記載の縦続接続型伝達系パラメータ推定方法であって、
   上記関数g_pは、
   

   

   
   、ただし、k=0,1,2,...,K、i=0,1,2,...,K、m=0,1,2,...,M、r=0,1,2,...,Mであり、K+M=P-1の関係を満たす
   ことを特徴とする縦続接続型伝達系パラメータ推定方法。
8. 過去の入出力に依存しない特性を持つ非線形系と有限長のインパルス応答特性を持つ線形系との縦続接続とみなせる伝達系について、当該非線形系のモデルパラメータおよび当該線形系のモデルパラメータを推定する縦続接続型伝達系パラメータ推定装置であって、
   Pをゼロ以上の整数として、時刻nにおける上記伝達系への入力信号x(n)に対し、P+1個の異なる関数g₀，g₁，g₂，…，g_Pを適用した結果と、当該各関数に対応する上記非線形系のモデルパラメータの推定値a^₀(n)，a^₁(n)，a^₂(n)，…，a^_P(n)との積和として、上記非線形系の出力推定値s^(n)を得る手段と、
   Lをゼロ以上の整数として、時刻n-Lから時刻nまでの各時刻に対応する上記非線形系の出力推定値s^(n-L)，…，s^(n)と、当該各出力推定値に対応する上記線形系のモデルパラメータの推定値h^₀(n)，h^₁(n)，h^₂(n)，…，h^_L(n)との積和として、上記線形系の出力推定値y^(n)を得る手段と、
   上記伝達系の出力信号y(n)から上記線形系の出力推定値y^(n)を引いた誤差信号e(n)を得る手段と、
   時刻n-Lから時刻nまでの各時刻に対応する各入力信号x(n-L)，…，x(n)を、P+1個の上記各関数g₀，g₁，g₂，…，g_Pに適用し、当該各関数毎に得られるL+1個の結果と、上記線形系のL+1個のモデルパラメータの推定値h^₀(n)，h^₁(n)，h^₂(n)，…，h^_L(n)との積和に相当する値を時間積和値として得、当該各関数毎のP+1個の上記時間積和値を要素として持つ非線形系入力評価ベクトルX(n)と、時刻n-Lから時刻nまでの各時刻に対応する各入力信号x(n-L)，…，x(n)を、P+1個の上記各関数g₀，g₁，g₂，…，g_Pに適用し、上記時刻n-Lから上記時刻n毎に得られるP+1個の結果と、上記非線形系のモデルパラメータの推定値a^₀(n)，a^₁(n)，a^₂(n)，…，a^_P(n)との積和に相当する値を関数積和値として得、当該時刻毎のL+1個の上記関数積和値を要素として持つ線形系入力評価ベクトルS(n)と、上記誤差信号e(n)、上記非線形系のモデル次数P、上記線形系のモデル次数Lを用い、Tをベクトルの転置として、上記非線形系のモデル次数Pおよび上記線形系のモデル次数Lに応じて、PがLよりも小さいほど更新量が大きくなるように正規化された上記非線形系のモデルパラメータの更新量で、上記非線形系入力評価ベクトルX(n)をスカラー倍したベクトルを上記非線形系のモデルパラメータベクトル[a^₀(n)，a^₁(n)，a^₂(n)，…，a^_P(n)]^Tに加えることで、時刻n+1における上記非線形系のモデルパラメータベクトル[a^₀(n+1)，a^₁(n+1)，a^₂(n+1)，…，a^_P(n+1)]^Tを推定し、上記非線形系のモデル次数Pおよび上記線形系のモデル次数Lに応じて、PがLよりも小さいほど更新量が小さくなるように正規化された上記線形系のモデルパラメータの更新量で、上記線形系入力評価ベクトルS(n)をスカラー倍したベクトルを上記線形系のモデルパラメータベクトル[h^₀(n)，h^₁(n)，h^₂(n)，…，h^_L(n)]^Tに加えることで、時刻n+1における上記線形系のモデルパラメータベクトル[h^₀(n+1)，h^₁(n+1)，h^₂(n+1)，…，h^_L(n+1)]^Tを推定する手段と
   を含む縦続接続型伝達系パラメータ推定装置。
9. コンピュータに、請求項１から請求項７のいずれかに記載の縦続接続型伝達系パラメータ推定方法の各処理を実行させるためのプログラム。
   

Images