WO2005066927A1

WO2005066927A1 - 多重音信号解析方法

Info

Publication number: WO2005066927A1
Application number: PCT/JP2005/000148
Authority: WO
Inventors: Shigeki Sagayama; Takuya Nishimoto; Hirokazu Kameoka
Original assignee: Toudai Tlo, Ltd.
Priority date: 2004-01-09
Filing date: 2005-01-07
Publication date: 2005-07-21

Abstract

多重音信号から発音数、基本周波数（F0）及び／あるいは周波数成分パワー比を検出する。複数の調波構造が混在したスペクトルのモデルを、単一の調波構造をモデル化した拘束つき混合正規分布モデルを混合することで構築する。本発明は二つの技術手段を含んでいる。一つは多重音信号の音源数推定であり、もう一つは各音源の基本周波数及び／あるいは周波数成分パワー比の推定である。拘束つき混合正規分布モデルのパラメータに関する推定と情報量基準に基づくアルゴリズムにより、発音数、それぞれの基本周波数及び／あるいは周波数成分パワー比を検出する。

Description

明細書

多重音信号解析方法

技術分野

[0001] 本発明は、多重音信号の解析法に係り、詳しくは、多重音信号から音源数を推定し、スペクトルを分離し、さらに複数の基本周波数 (F )を推定する問題を最適化問題とし

0

て定式ィ匕し、解を求める方法に関するものである。

背景技術

[0002] 音声や楽器音などの複数音源の音響信号が混在したものを多重音と呼ぶ。多重音の基本周波数推定は、音楽情報科学の分野にお!、て重要な要素の一つとして研究されている。効用としては、自動採譜、楽器音分離、演奏者の演奏分析、音響信号からの MIDI変換などが考えられる。

[0003] 従来より、音声の基本周波数推定法としてはケプストラム法や自己相関法などが用いられてきたが、音楽のように複数の基本周波数が混在する場合には適さな力つた。そこで、音楽情報科学の分野にぉ、て多重音または多声音楽を対象とした自動採譜の研究が 1990年頃力盛んになつたが、音源同士の基本周波数成分や倍音成分が互いに複雑に重なり合うことや、ミツシングファンダメンタル現象などが原因となり、音源数や各音源の基本周波数の推定は容易ではな力つた。

[0004] 柏野らは、和音構成音情報，和音遷移情報，楽音ごとの倍音比情報などの統計データによる仮説を統合したベイジアンネットワークを構成し，事後確率が最大となる仮説の組 (音名、音源名)、周波数成分を出力する手法を提案した。この手法では、観測パワースペクトルの極大点の集合の中から、各単音を構成する尤もらし、周波数成分を抽出することができるため、周波数解析における周波数分解能と同程度の分解能で基本周波数および倍音周波数を得ることができる。

[0005] また後藤は、単一音の倍音構造を確率分布としてモデルィ匕し、その確率分布モデルの混合分布の重みを推定する手法を提案した。これは、あらゆる基本周波数の単一音が同時発音をしていると見なして相対的に優勢な基本周波数を推定する手法であり、目的音の基本周波数を高精度に得ることができる。このように推定基本周波数値の分解能が高ければ高、ほど、音源の分離合成技術や会話音声などのように非音楽的な基本周波数推定の精度向上などが期待できる。

[0006] また、単一チャネル信号に対する多重 Fの推定の研究は、音楽情報科学以外にも、

0

雑音重畳下の音声認識、会議や討論などの状況を想定した同時複数音声認識、音声強調、韻律分析、音声符号ィ匕ゃ圧縮などに大きく貢献する。

[0007] Chazanらは、長い分析窓において時間に伴って連続的に変化する基本周波数に対して時間伸縮変換を最小二乗法により施すことで基本周波数を一定にしたのち、櫛形フィルタを用いて音声分離を行う手法を提案した。また Wuらは、フィルタバンク処理と、 Fダイナミクスを状態とした HMM (Hidden Markov Model)を用いた基本周波数

0

トラッキングによる雑音重畳化での多重 F推定手法を提案した。これらの手法は、精

0

度の高い基本周波数推定を実現し、良好な実験結果を得ているが、普遍性をもった音源数 (同時発話者数)判定基準にっ、ては議論されて!、な、。

特許文献 1：特許第 3413634号

非特許文献 1 :柏野邦夫，中臺一博,木下智義，田中英彦: "音楽情景分析の処理モデル OPTIMAにおける単音の認識,"電子情報通信学会論文誌， D-II,

Vol.J79-D-II,No.l l, pp.1751- 1761, 1996

非特許文献 2 :柏野邦夫，木下智義，中臺一博，田中英彦： "音楽情景分析の処理モデル OPTIMAにおける和音の認識,"電子情報通信学会論文誌，

Vol.J79-D-II,No.l l, pp.1762- 1770, 1996

非特許文献 3 :後藤真孝:"音楽音響信号を対象としたメロディーとベースの音高推定, "電子情報通信学会論文誌， D-II, Vol.J84-D-II,No.l, pp.12-22, 2001

非特許文献 4 : M. Goto: "A Predominent-F0 Estimation Methodfor Real-world Musical Audio Signals: MAP Estimation for Incorporating Prior Knowledgeabout FOs and Tone Models, "Proceedings of CRAC— 2001},2001

非特許文献 5 : M. Goto: "A Predominent-FO Estimation Methodfor CD Recordings: MAP Estimation Using EM Algorithm for Adaptive Tone Models , " Proceedings of the 2001 IEEE International Conference on Acoustics, Speech and Signal Processing, ICASSP〜2001, pp.V- 3365- 3368, 2001 特許文献 6 : D.Chazan,Y.Stettiner and D.Malah, "Optimal Multi-pitch Estimation Using the EM Algorithm for Co- channelSpeech Separation, 'Proc. ICASSP93}, Vol.2, pp.728— 731,1993

非特許文献 7 : M.Wu, D.Wangand G.J. Brown, " AMulti-pitch Tracking Algorithm for Noisy Speech,"lCASSP2002},Vol.1 , pp.369- 372, 2002

発明の開示

発明が解決しょうとする課題

[0008] 本発明は、多重音信号から複数の基本周波数 (F

0 )を抽出することを目的とするものである。本発明のさらなる目的は、多重音信号の基本周波数抽出における音源数（本明細書において、音源数とは、同時発音する基本周波数の数をさす。）を抽出することにある。本発明のさらに他の目的は、多重音信号から内在する倍音構造ごとにスペクトルを分離し、周波数成分パワー比の抽出を行うことである。

課題を解決するための手段

[0009] 本発明は、多重音信号の観測スペクトルを、複数の調波構造モデルからなる多重調波構造モデルで推定することで多重音信号を解析する方法である。各調波構造モデルは、基本周波数成分に対応する一つの単峰形分布の代表値である基本周波数推定値と、該基本周波数推定値によって決定される他の単峰形分布の代表値とを有する。単峰形分布としては数々の分布が知られている力一つの好ましい態様では、単峰形分布は正規分布 (対数正規分布を含む)である。分布の代表値としては、平均、中央値、最頻値が例示されるが、一つの好ましい態様では、分布の代表値は平均である。後述する実施形態では、多重音信号の観測スペクトルを、拘束つき混合正規分布によりモデル化した調波構造モデルを複数混合した多重調波構造モデルで推定する多重音信号の解析方法について説明するが、拘束つき混合正規分布によりモデル化した調波構造モデルは、基本周波数成分に対応する一つの正規分布の平均である基本周波数推定値と、該基本周波数推定値によって決定される他の正規分布の平均とを有している。

[0010] 多重調波構造モデルにおけるモデルパラメータは、調波構造モデルの各単峰形分布の代表値、重み、分散から成る。代表値パラメータは、基本周波数推定値を構成する代表値を含む調波構造モデルの各単峰形分布の各代表値により構成されるが、代表値パラメータにおいては、基本周波数推定値のみが自由パラメータであり、他の代表値は基本周波数推定値によって拘束されるパラメータである。重みパラメータは、調波構造モデルごとの重み及び調波構造モデルにおける単峰形分布の重みの積として表される。分散パラメータは、一つの態様では、既知パラメータとしてモデルに与えられる。

[0011] 本発明は二つの技術手段を含んでいる。一つは多重音信号の音源数推定であり、もう一つは各音源の基本周波数及び Zあるいは周波数成分パワー比の推定 (調波成分エネルギー検出）である。二つの技術手段は互いに独立して成立するものである力二つの技術手段を組み合わせることで、多重音信号から音源数を推定し、スぺクトルを分離し、さらに複数の基本周波数 (F )を推定することができる。さらに、後者の

0

技術手段は、基本周波数の検出と、周波数成分パワー比の検出という二つの技術手段を含んでヽる。基本周波数検出と周波数成分パワー比検出は検出ステップにお!/ヽては共通するが、それぞれ独立の技術手段として採用することが可能である。

[0012] 多重音信号の音源数推定は、モデルパラメータを推定するにあたり、所与の一つ又は複数の調波構造モデルに基づ、てモデル数を減少あるいは Zおよび増力！]させながら情報量基準を用いて適切な複数の調波構造モデルを選択し、選択された調波構造モデルの数を推定音源数とするものである。モデル数の減少には、モデルの削除やモデルの統合が含まれる。モデル数の増加には、モデルの分裂や新たなモデルの生成が含まれる。一つの態様では、多重音信号の音源数推定は、推定される音源数より多い数の調波構造モデルを用意し、尤度への寄与度が少ないと判断される調波構造モデルを削除しながら情報量基準を用いて適切な複数の調波構造モデルを選択する。尤度への寄与度が少ないと判断される調波構造モデルは、調波構造モデルの調波構造ごとの重みが小さいモデルである。より具体的には、すべての調波構造モデルの中で重みが最小のモデル、あるいは、隣接する基本周波数推定値が所定の閾値より近接した場合には、重みが小さい方のモデル、のいずれかに該当するモデルが削除される。この場合、調波構造モデルにおける各単峰形分布の重みを固定しておいても良い。後述する実施例では、情報量基準として AICを採用し、 AIC が最小となるときの調波構造モデル数を推定音源数として、る。調波構造モデル数を減少させながら適切なモデルを選択することにつ、て説明したが、調波構造モデル数を増加させながら適切なモデルを選択するものでもよい。例えば、ある調波構造モデルの代表値が他の代表値力孤立しているものは分裂させた方が良い場合がある。一つの基準としては、「特定のモデルの代表値力それと最も近くにある他のモデルの代表値との距離がある一定閾値より大きい場合に、その特定のモデルを 2つに分裂する」ことが例示される。この場合、もし分裂したことのメリットがな力つた場合には 2つの代表値は再び引きつけ合い、また統合されるし、 2つに分裂したことで多峰分布がうまく表現できた場合には適切な距離を保って収束するのではと推測される。

[0013] 基本周波数あるいは Zおよび周波数成分パワー比の検出は、推定音源数の基となつた各調波構造モデルにぉヽて、基本周波数推定値を構成する代表値を少なくとも整数倍および Zあるいは整数分の 1倍の値に順次置き換えながら少なくとも調波構造モデルの各単峰形分布の重みを推定するにあたり、情報量基準を用、て適切な置き換え値及び重みを選択し、選択された置き換え値あるいは Zおよび重みを、基本周波数あるいは Zおよび周波数成分パワー比の推定値とするものである。一つの態様では、基本周波数推定値の置き換え値は、基本周波数推定値の少なくとも 1倍

、 2倍、 3倍、 2分の 1倍、 3分の 1倍、 3分の 2倍、 2分の 3倍からなる群から選択された 2つ以上の値である。一つの態様では、パラメータの最尤推定量は、基本周波数推定値 (推定代表値)を置き換え値に固定して単峰形分布の重みのみを更新しながら求める。基本周波数は、置き換えられた基本周波数推定値（1倍を含む)から求め、周波数成分パワー比は、単峰形分布の重み力求める。後述する実施例では、情報量基準として AICを採用し、 AICが最小となるときの置き換えられた基本周波数推定値を推定基本周波数、分布の重みを推定周波数成分パワー比とする。

[0014] 本発明は、推定手法と情報量基準を用いるものであり、後述する実施例では、一つの好ま、態様として、 EMアルゴリズムを用いたモデル/スペクトル間の K-L情報量最小化によるモデルパラメータ推定 (最尤推定法に相当する）と両分布間の二乗誤差に基づく対数尤度を用いて計算される AICとの組み合せについて説明する。しかしながら、多重調波構造モデルにおけるモデルパラメータの推定手法としては、最尤推定法の他に、最大事後確率推定法 (MAP推定)、変分ベイズ法が含まれる。また、一つの態様では、モデルパラメータの推定は、観測スペクトルと多重調波構造モデル間の K-L情報量を評価尺度として行、、モデルパラメータの最尤推定量を K-L情報量を最小化するパラメータで近似する。本発明が採用する情報量基準についても AIC (赤池情報量基準）に限定されるものではなぐ ABIC (赤池ベイズ型情報量規準 )、 TIC (竹内情報量基準)、 BIC (ベイズ情報量基準)、 MDL (最小既述長）、 GIC ( 一般化情報量基準、正則化法 (罰則つき最尤法)、 EIC (ブートストラップ情報量基準 )、 PIC (予測情報量基準）、クロスノくリデーシヨン、 FPE (最終予測誤差）、マローの C 基準、ハナンクインの基準、さらに、これらの近似的な情報量基準あるいは等価の

P

情報量基準を含み、例えば、 MAP推定と BICあるいは MDLの組み合わせや変分ベイズ法 (モデル選択を同一の枠組として含む）を採用することができる。また、モデルパラメータの推定に用いる推定アルゴリズムとしては EMアルゴリズム以外の既知の手法、例えば、最急降下法を EMアルゴリズムに代えて、あるいは、 EMァルゴリズムと組み合わせて採用することができる。また、本発明は、本発明に係る信号解析方法をコンピュータに実行させるためのコンピュータプログラム、あるいは、このようなコンピュータプログラムを記録させた記録媒体としても提供される。

発明の効果

[0015] 本発明によれば、音源数が未知の多重音信号力音源数 (発音数)を検出すること力 Sできる。また、多重音信号力も基本周波数 (音の高さ)を精密な値として検出することができる。また、多重音信号スペクトル力各調波構造を分離し、周波数成分パヮ一比を検出することができる。

発明を実施するための最良の形態

[0016] [A]拘束つき混合正規分布モデルの最尤推定

[A - 1]多重調波構造モデル

短時間スペクトルの解析では、解析区間に窓関数を掛けることが一般的である。そのため、周波数が一定の単一正弦波の信号であっても、線スペクトルではなぐ左右に広がりをもつスペクトルが観測される。これは、窓関数のフーリエ変換と線スペクトルとの畳み込みを行うことに相当するためである。さらに分析窓区間内で周波数が連続的に変化する場合、それに応じた広がりをもつスペクトルが観測されることになる。これらに起因する基本周波数成分や調波成分の広がりにより、異なる信号同士の周波数成分が重なり合い、近接する周波数成分の分離や正確な基本周波数 (F )あるい

0 は高調波周波数の検出が困難となる。

[0017] このように広がって観測される周波数成分を各周波数の出現頻度分布あるいは確率分布と見なし、その分布を正規分布により近似することで、単一の調波構造を有するスペクトルを複数の正規分布の混合分布としてモデルィ匕する。図 2に示すように、スベクトルの拡がり形状を正規分布で近似することで、周波数値を正規分布の平均推定、周波数成分のエネルギーを混合正規分布の重み推定に対応させることができる。調和性の保持のため、基本周波数成分に対応する 1つの正規分布の平均 (これを以後基本周波数推定値と呼ぶ)のみが自由度をもち、その位置に応じて残りのすべての正規分布の平均の位置は決定される。単一の調波構造をこのような拘束つきの混合正規分布によりモデルィ匕したものを本明細書にぉ、て「調波構造モデル」と呼ぶことにする（図 3参照)。基本周波数推定値を； z と置けば、調波構造モデル kの各 k

平均 kは、対数周波数領域において、

[数 1] fc = i fc， fc + log 2，'，'， _fc +log n，- ' ' , _fc + log iV_fc} (1) のようなベクトルで表される。ただし、 nは n次高調波成分に対応する正規分布のラベルを、 Nは調波構造モデルごとの Nyquist周波数まで取り得る正規分布の平均の数を k

表す。

[0018] 複数の調波構造が重なり合うスペクトル（図 1参照）を、調波構造モデルをさらに混合することによりモデルィ匕し、これを多重調波構造モデル P (X)と呼ぶことにする。ただし

Θ

、 Xは対数周波数とする。 P (X)のモデルパラメータ { Θ }は、

Θ

[数 2]

W = (M^ Wk ^ I & = 1, · · -，Κ} (2) であり、 w^k={w^k , · · · , w^k , - - -w^k }は調波構造モデル kの各正規分布の重み、 σは分散（σは固定であり、実験においては予め与えられる定数とする)、 Κは混合された調波構造モデルの数をさす。

[Α— 2] ΕΜアルゴリズムによるモデルパラメータ推定

正規ィ匕した観測スペクトル χ)と上記した多重調波構造モデル Ρ (X)との

Θ

Kullback- Leibler情報量 D( θ )は以下となる。

[数 3]

^DW (3)

D( θ )を最小化するモデルパラメータを求めることと、モデル Ρ (X)の平均対数尤度

Θ

( χ)を対数周波数 Xの出現頻度を表す確率分布と解釈した場合、対数尤度の Xに関する期待値)を最大化するモデルパラメータを求めることは等価である。そこで、 Xについて、モデルパラメータ Θを Θ (バー)に更新したときのモデルの平均対数尤度の差は画 (^) logP^(x) -/( ) log ,(2:) = /( 。g^ (4)

ΡΘ {Χ)

となる。 Dempsterらは、式 (4)において x)を確率密度分布関数とし、平均対数尤度を最大にするために ΕΜアルゴリズムを定式ィ匕した。 Xがどの正規分布によって生成されたかは一意に決定できないため、これを直接最大化することはできない。そこで両辺に対し、 Xがどの正規分布力生成されたかを表す P (n, k|x)についての期待値を

Θ

とることで Q関数と呼ぶ以下のような評価関数

[数 5]

Q{e, 6) {x) log P_d (x,_n,k)dx (5)

を導出することができ、

[数 6]

が成立するため、 Q関数を最大化する Θ (バー)を求めることで、 Xに関するモデルの平均対数尤度を単調に増カロさせることができる。 g(x|x

0， σ ²)を平均 X、分散 σの正規

0

分布とすると、 P (n, k|x)は、

Θ

[数 7]

と書け、 logp χ, _n, k)は各正規分布が与える対数尤度なので、

Θ (バー )（

[数 8]

である。

[0020] 以上より、初期設定 (ステップ 0)を経て、以下のような Εステップと Μステップによる反復計算の収束性は保証され、モデルの K-L情報量を局所最小化するパラメータ , wを得ることができる。

[0021] (1)ステップ 0:

(初期設定）

モデルパラメータ μ , wの初期値を与える。 (2) E-ステップ：

式 (5)により Q( 0 , 0 (バー))を計算する。

(3) M-ステップ：

Q( θ , Θ (バー》を最大化する Θ (バー)を計算する。

[数 9]

Θ = &rgma.x Q{6, S) (11)

Θ

Θを更新後、 Ε-ステップに戻る。

尚、正規分布に代えて対数正規分布を用いる場合であっても、実質的に同様のステツプを適用することができる。また、本発明に用いるスペクトル解析には、離散的フ一リエ変換 (FFTを含む）によるスペクトル計算、ウェーブレット変換によるスペクトル計算が例示される。

[0022] [Α— 3] Clusteringとしての解釈

この拘束つき混合正規分布モデルの最尤推定は、スペクトル密度分布を架空の微小エネルギーの度数密度分布と捉えた場合、微小エネルギーを Clusteringにより各音へ分類する問題であると解釈することもできる。

[0023] 正規分布の平均 μ k+lognをクラスタ n、 kの中心と考えれば、事後確率 P (n, k|x)を xに

Θ

位置する微小エネルギーがクラスタ n、 kに帰属する確率、対数尤度 log P (n, k, x)を

Θ

クラスタ kの中心と χに位置する微小エネルギーとの距離を表す関数と見なすことができる。 Clusteringの評価関数は一般に、各微小エネルギーの位置と帰属するクラスタ中心との距離の自乗の和で表されるため、上記の観点により Q関数と同一と見ることができる。ただし、 Q関数では事後確率 P (n, k|x)と対数尤度 log P (n, k, x)は同一

Θ Θ

の確率分布に基づいて計算されるのに対し (いずれも g(x k+logn, σ ²)を用いて計算される）、 Clusteringにおいては必ずしもクラスタ帰属確率とクラスタ中心との距離関数は同一の関数 (正規分布など）を対応させる必要はない。例えば、微小エネルギーを最近傍のクラスタ中心のクラスタにすべて帰属させ、距離関数をユークリッド距離の自乗とした場合、帰属確率と距離関数は全く別の関数で表現される。このとき評価関数の最大化問題は、 k-meansアルゴリズムとして定式化できる。このように、多重基本周波数 (F )推定問題として定式化される前述した拘束つき混合正規分布の最尤推

0

定を包含する、 Clusteringの考え方による多重基本周波数 (F )推定手法を「

0

Harmonic Clusteringと呼んでいる。クラスタリングの考えを採用した場合には、前述の「調波構造モデル」は「倍音クラスタ群 (倍音構造としての拘束をもつ複数のクラスタ）」、「平均」は「クラスタ重心」、「基本周波数推定値」は「基本クラスタ重心」と言い換えることができる。

[0024] [B]多重基本周波数 (F )推定アルゴリズム

0

本発明に係る多重基本周波数 (F )推定アルゴリズムは二つの処理プロセスにより構

0

成される。音源数 (含む、同時発話者)推定検出プロセス、及び、基本周波数 (F )お

0 よび周波数成分パワー比検出プロセスである。これらのプロセスは情報量基準に基づいて検出される。上記までは K-L情報量最小化のステップについて説明した。今、モデルとスペクトル間の偶然誤差が正規分布型の確率分布に従うと仮定した場合、偶然誤差の対数尤度を最大化する最尤パラメータは本来最小二乗誤差パラメータに相当する力上記のように両分布間の K-L情報量を最小化するパラメータであれば E Mアルゴリズムにより効果的に推定できるので、近似的にこれを疑似最尤パラメータと考え、両分布間の二乗誤差総和を計算して疑似最大対数尤度を求めることができる。しかし、モデルの比較基準として最大対数尤度は必ずしも適切ではない。そこで、複数のモデルが候補にあるとき、適切なモデルを選択する基準として情報量基準を用いる。情報量基準としては、 AIC (赤池情報量基準）、 ABIC (赤池ベイズ型情報量規準)、 TIC (竹内情報量基準)、 BIC (ベイズ情報量基準)、 MDL (最小既述長)が例示されるが、一つの好ましい態様では、情報量基準として AICが用いられる。

[0025] [B-1]AICによるモデルの選択基準

1つのモデルに対する最尤パラメータを求めることができ、自由パラメータ数に応じてモデルの候補が複数個あるとき、その中から適切なモデルを選択する規準が必要となる。そこで、赤池によって提唱された AIC(Akaike Information Criterion,赤池情報量規準)を導入する。 AICは、

[数 10] AIC = -2 x (モデルの最大対数尤度）

+2 x (モデルの自由パラメータ数） ²) で与えられ、適切な自由パラメータ数のモデルを選択する問題において有効であることが知られている。

[0026] [B-2]音源数（同時発話者数)検出プロセス

EMアルゴリズムにより得られるモデルパラメータの最小 K-L情報量は初期値に依存し、しばしば誤った局所解に陥る場合がある。そこで、基本周波数推定値 kの誤つた局所解への収束を回避するため、予想される音源数（同時発話者数)より多めの数の調波構造モデルを満遍なく初期配置しておくことで目的とする解が得られる可能性は高くなるはずである。ただし、このように初期配置された調波構造モデルの数が音源数（同時発話者数)より多ぐかつすベての目的解が得られているならば、多重調波構造モデルは観測スペクトルに対して明らかに過適応を起こしていると言える。もし、周波数成分の分布の形状が正規分布で十分近似可能であれば、調波構造モデルは音源数（同時発話者数)と同数あれば十分なはずである。そこで、不必要な調波構造モデルを順次削減して、き、 AICが最小となる調波構造モデル数を判定することで音源数（同時発話者数)を推定する。具体的な処理手順を以下に示す。

[0027] (1)任意の周波数区間に基本周波数推定値 1, · · · , k}の初期値を設定する。

(2) [A-2]で述べた EMアルゴリズムにより最小 K-L情報量を求める。ただし、ここでは正規分布の重み w ^kに関して

[数 11]

という拘束を与える。これは、正規分布ごとではなく調波構造モデルごとの重みを規準として、優先的に削減すべき調波構造モデルを決定するためである (後述)。この場合、 M-ステップにおける , w^kの更新値は式 (5)のそれぞれのパラメータに関する

k

偏微分を 0と置くことで得られ、以下で与えられる。

Pn (^)f(^) dx

(3)最小 KL情報量パラメータにおけるモデルとスペクトル間の二乗誤差総和 (積分）により求めた最大対数尤度を用いて AICを算出する。調波構造モデルごとに 2つの自由パラメータ μ , w^kがあるので、自由パラメータ総数は 2 X Kである。 AICが最小と k

なるときの調波構造モデル数を推定同時発話者数とする。

(4)以下の、ずれかを満たす調波構造モデルを削減し、残った調波構造モデル数を K (バー）とする。 K= Κ (バー）として（2)に戻る。

[0028] (4-1)すべての調波構造モデルの中で重み w^kが最小のもの。これは、重み最小の調波構造モデルは全体のモデルが与える期待対数尤度に及ぼす重要度が低いと見なせるためである。

(4-2)隣接する 2つの基本周波数推定値がある一定閾値より近接した場合、 w^kが小さい方。これは、 1つの極値に 2つの基本周波数推定値が収束していると考えられるためである。

すなわち、情報量基準に基づき適切な調波構造モデル数となるまで重要度 (調波構造ごとの重み w^k)の低、ものから順に削除して！、く。

[0029] このプロセスを図 4のスペクトルに対して実際に行った例を図 5に示す。図 5の下図における「+」は基本周波数推定値の反復計算ごとの更新値、破線が（ 2 )において収束判定によりモデルパラメータ/ z , wが収束したと見なされた時点を表す。上図の折れ線グラフが各時点での AICの値を表す。調波構造モデル数が 3のときに AICは最小値をとるため、この場合推定同時発音数は 3となる。

[0030] [B— 3]基本周波数 (F )および周波数成分パワー比検出プロセス前述したプロセスにおいては、基本周波数推定値が真の基本周波数 (F )とその k 0 整数倍あるいは整数分の 1倍の値のときも K-L情報量を極小にすると考えられるため、得られる基本周波数推定値は必ずしも真の基本周波数 (F )であるとは限らない。

0

そこで、前述のプロセスで得られる基本周波数推定値の最小 K-L情報量推定値 k

1S 真の基本周波数 (F )あるいはその整数倍か整数分の 1倍の値であることを前提と

0

し、 μ をそれぞれの値に順次置き換えながら何らかの規準に基づいて真の基本周 k

波数 (F )を検出する。

0

[0031] ここで、前節で与えた重み w ^kに関する拘束を外し、すべての正規分布の重みに関して自由度を与えることにする。重みの総和は 1なので、厳密には、自由度は 1つを除く残りすベての重みに対して与えられる。従って、重み w ^kの推定値は近似された周波数成分パワー比 (調波成分間の相対的な強度)を表すことになる。もし、置き換えたが真の基本周波数 (F )より小さい場合、実際に存在する調波成分に対応する k 0

ガウス成分以外の成分が実際に存在してもしなくてもモデルとスペクトル間の二乗誤差総和に変化はほとんど及ぼさないと考えられるため、過適応を起こしていると言える。例えば、 μ が真の基本周波数 (F )の 1/2に対応した場合に重み w ^kに関して Ε k 0 n

Mアルゴリズムを行えば、偶数次の調波成分に比べて奇数次の調波成分が極端に小さい単一音のモデルとして表現されるはずである (図 6)。また、 μ が真の基本周波 k

数 (F )の整数倍の場合には、真の基本周波数 (F )と比較して二乗誤差の総和 (対

0 0

数尤度）は十分に小さいと考えられる。この観点から、前節のプロセス同様、 AICに基づいて真の基本周波数 (F )を検出することができると考えられる。前節のプロセス

0

にお、て残った調波構造モデルすべてにつ!、て以下の手順を行、、基本周波数 (F )および周波数成分パワー比の検出を行う。尚、次に述べる手順では、後述する実

0

験を行う上での便宜のため、基本周波数推定値を整数倍したものの中カゝら基本周波数 (F )を探索しているが、基本周波数推定値を整数分の 1倍したものの中からも基

0

本周波数 (F )を探索することが望ましい。

0

[0032] (1)調波構造モデル kにおける基本周波数推定値を +logtに置き換える。ただし、 t k

は初期値 1の整数とする。このとき、上限がナイキスト周波数の対数である範囲内にとりうる正規分布の数を N する。 (2) EMアルゴリズムにより最小 K-L情報量パラメータを求める。ここでは、更新すべきパラメータは各正規分布の重み w ^kだけである。 M-ステップにおける w ^kの更新値は以下で与えられる。

[数 13] =ゾ ( dx (16)

(3)最小 Κ-L情報量パラメータにおけるモデルとスペクトル間の二乗誤差総和 (積分）により求めた最大対数尤度を用いて AICを算出する。このとき、自由パラメータ総数は N ¹である。 tを 1増やし、 ( 1 )に戻る。 AICが最小となるときの +logtを推定基本周 k k

波数 (F )とする。また、この時の、重み w ^kの最尤推定量が、推定周波数成分パワー

0 n

比に相当する。

[0033] これらを統合した基本周波数抽出アルゴリズムとしての性能を、モノラルの音楽音響信号を対象とした実験により評価する。本発明に係る信号解析はコンピュータを主要構成とする信号解析装置によって行われ、該信号解析装置は、信号取り込み部、信号解析部 (CPU等の処理装置)、解析結果等を表示する表示部 (表示画面)、記憶部 (メモリ及び外部記憶装置を含む）、入力部（マウス、キーボード等）、コンピュータを動作させる制御プログラム等を備えて、る。

実験例 1

[0034] 2段階の処理行程を統合した基本周波数抽出アルゴリズムの動作実験を、モノラル音楽音響信号に対して行った。スペクトル解析は、サンプリング周波数を 44.1kHz,フレーム長を 25ms,フレームシフトを 10msとし、 Hamming窓を窓関数として FFT (高速フ一リエ変換)を行い、短時間スぺ外ル系列を得た。実験に用いた対象曲の曲名、演奏方式、演奏者は図 7の通りである。

[0035] 性能評価の目安として以下のヒューリスティックな方法により音名正解率を与えた。抽出した基本周波数を最も近、音名の基本周波数に割り当てる。実験データの音響信号を 1小節ごとに分割し、各小節のフレーム数を得る。楽譜に記譜されている音符に応じて各音の音長フレーム数 (四分音符の場合は 1小節のフレーム数の 1/4)を与え、音長フレーム数の 1/2以上のフレーム数分だけ正解音名に相当する基本周波数を抽出できた場合に、 1つの音に対して正解とした。また、同一のピッチクラスで、オタターブ位置の異なる複数の音が同時発音している場合には、オクターブ位置の低い方の基本周波数を抽出できていれば、正解とした。正解した数を以上の基準に従って目視により数え、総音符数を分母として正解率を計算した。

[0036] ヴァイオリン三重奏による Kanonに対する基本周波数の抽出結果の一部をそれぞれ楽譜とともに図 8に示す。 X印がフレームごとに抽出した基本周波数値を、破線が正解音名に相当する基本周波数を表す。音名正解率は、 92.7%であった。図 8の通り、全体として音源数およびオクターブ位置に関して妥当な推定結果が得られていることが分かる。

実験例 2

[0037] 上述した基本周波数 (F )検出アルゴリズムの性能を確認するため、基本周波数 (F )

0 0 検出方法としてよく知られる Cepstrum法と比較実験を行った。 ATR音声データベースより音声データとハンドラベルによる基本周波数（F )パターンの referenceデータを用

0

いた。すべての音声信号はサンプリング周波数 12kHzでディジタルィ匕され、フレーム長 64ms、フレームシフト 10msのもとで Hamming窓をかけて周波数解析 (FFT)を行!、、スペクトル系列を得た。同時発話者数検出プロセスにおいて、初期調波構造モデル数は 4とし、基本周波数推定値を配置する周波数範囲は 70Hzから 140Hzとした。また、すべての正規分布の分散の値は 0. 45とした。 fymおよび myiから始まる音声ファイル名はそれぞれ女性話者と男性話者による音声信号データをさす。評価基準として、検出された基本周波数 (F )が基本周波数 (F )パターンの referenceデータの値から

0 0

5%以上外れた場合は、 gross errorと見なした。

[0038] 話者一人による音声信号に対する実験について説明する。提案手法が多重基本周波数 (F )についてだけではなく単一基本周波数 (F )についても高い性能で推定で

0 0

きることを確認するため、話者一人による単一チャネル音声信号に対して動作実験を行い、単一基本周波数 (F )推定手法として広く知られる Cepstrum法と推定正解率の

0

比較を行った。推定正解率 (Accuracy)は、総フレーム数に対する gross error以外のフレーム数の割合とした。 [0039] 女性話者および男性話者それぞれの音声データにおける実験結果を Cepstrum法の結果と併せて図 10に示す。また、基本周波数お )検出結果の例を図 9上図に示し、

0

対応する reference基本周波数 (F )パターンを図 9下図に示す。結果より、推定正解

0

率 92. 4%— 99. 0%を得た。また、すべての音声信号に対して Cepstrum法に比べて推定正解率が高ぐ多重基本周波数 (F )だけではなく単一基本周波数 (F )推定手法

0 0 としても高、性能であることが確認できた。

[0040] 次に、話者二人による同時発話音声信号に対する実験について説明する。話者二人による単一チャネル音声信号に対して動作実験を行、、同様に Cepstrum法との比較を行った。 C印 strum法は複数話者による発話には原理的には適用できな、ため、厳密には客観評価の比較対象とはならないが、提案手法の客観評価のための参考基準としては十分であると考えた。

[0041] 2つの音声データの信号波形を人工的に加算したものを同時発話音声データとし、 SSR(signal-to-signal ratio)は OdBとした。 Cepstrum法による基本周波数（F )検出は、

0 低ケフレンシ一領域と高ケフレンシ一領域を閾値により区分し、高ケフレンシ一における 2つのローカルピークを抽出することで行った。推定正解率は、同時発話されていると見なされるフレームを referenceFパターンから判断し、同時発話時のフレーム

0

総数に対する gross error以外のフレーム数の割合とした。

[0042] Cepstrum法の推定正解率を図 13、提案手法の推定正解率を図 14に示す。また、提案手法の基本周波数 (F )検出結果の例を図 11, 12の上図に示し、それぞれに対応

0

する referenceFパターンを図 11, 12の下図に示す。 Cepstrum法では、推定正解率

0

力 S41. 0%— 71. 7%程度であつたのに対し、提案手法では、 72. 6%— 92. 8%であった。同時発話者数を事前に与えなくても AICにより基本周波数 (F )および話者数を高精

0

度に推定することができ、情報量規準が多重基本周波数 (F )推定に十分有効である

0

ことが確認できた。

[0043] [C]拘束つき混合正規分布モデルの MAP推定による同時発話音声の基本周波数追跡及び調波構造分離

[C—1]拘束つき混合正規分布モデルの定式ィ匕

短時間周波数解析では、一般に解析区間に窓関数を掛けるため、左右に広力 ^をもつスペクトルが観測される。窓関数として正規分布窓を用いれば、スペクトルの広がりの形状が理論的に正規分布の形状となるので、基本周波数成分に対応する正規分布の平均だけが自由度をもつ拘束つき正規混合分布により単一音の調波構造をモデルィ匕できる。 k番目の調波構造モデルの各平均は、

[数 14]

· · · , ημ/c , ' · · , と書ける。ただし、 ηは η次高調波成分に対応する正規分布のラベルを、 Νは正規分

k

布の数を表す。

K個の音の調波構造が重なり合うスペクトルを、単一の調波構造モデルを K個混合することによりモデル化し、モデルパラメータを、

[数 15]

とする。 w^k , σ ^kは n次成分の重み、分散を表す。スペクトル分布を正規化して確率変数 (周波数） ωの確率分布 f ( co )とみなせば、 Θの事後確率を最大化する Θは、以下の式で表される。

[数 16]

Θ (17)

ρ ( θ ) « θの事前確率を表す。 ΕΜアルゴリズムにより以下の Q関数を最大化する Θ (ハット)を 0の更新値として逐次的に計算することで局所最適解を得ることができる。

[数 17] Q (

[0044] [C-2]同時発話音声の基本周波数追跡アルゴリズム

話者数および各基本周波数を検出する処理 (検出処理)と直前フレームにおいて検出された基本周波数に基づき基本周波数を追跡する処理 (追跡処理)のヽずれか一方をフレーム毎に実行し、逐次的に複数の基本周波数を同時検出していく。発話開始時、フレーズ境界や新たな話者の音声介入時などにおいては、話者数とそれぞれの基本周波数を検出する必要がある。初期フレームは「検出処理」を実行し、以降のフレームでは、直前フレームでの多重音モデルと f ( co )の KL(Kullback-Le¾ler)情報量が一定閾値以下の場合は「追跡処理」を、閾値より大き!、場合は新たなフレーズの開始直後あるいは新たな話者による音声の介入直後と見なして改めて「検出処理」を行う。検出処理は、話者数推定ステップと基本周波数検出ステップカゝら成る。

[0045] [C 3]話者数推定ステップ

EMアルゴリズムにより得られるモデルパラメータの収束値は初期値に依存し、しばしば誤った局所解に陥る。そこで、 μ が目的解へ局所収束する可能性を高くするため

k

、予想される発音数よりも多めの調波構造モデルを満遍なく初期配置する。ここで、 AICを導入し、適切な自由パラメータ数を推定する。すなわち、不必要な調波構造モデルから削減していき、 AICが最小となるときの数を推定話者数とする。具体的な手順を以下に示す。

(1)限定した周波数帯域内に基本平均を K個配置する。

(2) EMアルゴリズムにより Θの最尤推定値を求める（事前分布を一様分布とする）。ここで、正規分布の重みは kのみに依存する調波構造モデルごとの重みパラメータ w kとする。尚、調波成分の強度比を事前にモデルに与えることも可能である。式（18) を最大化する， w^k, a ^kの更新値は偏微分を 0と置くことで以下のように求められ

k n

る。

[数 18]

— oo

J一 oo P(⁷ Ι^ω,の ― ημ„)²άω

ΧΠο Ρ(^η， ^,の/ (w)d j

(3) AICを算出する。 AICが最小値をとるときの調波構造モデルの数 K (ハット）を推定発話者数とする。

(4) w^kが最小の (尤度への関与が最も低ぐ不必要とみなせる)調波構造モデルを削除する。 σ ^kを大きめの値に置き換えて、ステップ（2)に戻る。尚、分散 σ ^kの更新は、分散の推定値を得るためではなぐ大きい初期値を与えることで kの目的解への収束を促進させるために行う。

[C 4]基本周波数検出ステップ

前ステップにおいて求められるの局所最適解は、真の基本周波数だけではなくそ k

の整数倍あるいは整数分の 1倍のいずれかに対応する可能性がある。ここでは、各調波成分の強度を手がかりとして真の基本周波数を検出する。 μ

kを整数倍/整数分の 1倍に置き換えながら、その都度正規分布ごとの重み w^kの最尤推定値から調波成分の強度比を推定する。 AICに基づき、適切なの位置を推定する。前ステップ k

にお、て残った調波構造モデルすべてにつ!、て以下のステップを行う。

(1)調波構造モデルの 1次成分の平均を t に置き換える。但し、 tを初期値 1の自 k

然数とする。限定した周波数帯域内までとり得る正規分布の数を Ntkとする。

( 2) EMアルゴリズムにより w^k , σ ^kの最尤推定値を求める。 Mステップにおける更新値は、式（22)、式（21)となる。

[数 19] p(n, k\ , θ) άω

(3)自由パラメータ数をとして AICを算出する。 tを 1増やし、ステップ（1)に戻る。

k

AICが最小となる t が推定基本周波数となる。また、最終的な w^kの最尤推定値が

k n

各音の周波数成分パワー比 (調波成分強度比)の推定値となる。

[C 5]追跡処理

1つのフレーム区間では、ある時点の基本周波数と直前の基本周波数との間には強い依存関係があるはずである。そこで、直前フレームでの基本周波数の検出結果を μ の事前分布に反映させ、最大事後確率 (MAP)推定によりをフレーム毎に更 k k

新 (追跡)する。 IX

kの（直前フレームでの

kの推定値に基づく）予測値を（ハット） k とし、 μ の事前分布を、 μ (ハット）を平均、 Vを分散とした正規分布とすれば、式（1 k k

8)より EMアルゴリズムの Mステップにおける μ の更新値は、

[数 20]

となる。また、重み w^kと分散 a ^kの更新はそれぞれ式 (22)、式 (21)を用いる。この追跡処理が連続で 3フレーム以上続く場合、予測値（ハット）は、過去の直前の 2フ k

レームにおける; z の推定、を結ぶ直線の延長上の値とし、 μ (ハット） = 2 k k k k β —β kと定める。それ以外の場合は、直前フレームの推定値を予測値とする。

[C-6]調波成分の周波数成分パワー比パラメータの導入

音 kの n次調波成分と基本周波数成分との強度比を r^k (r^k = 1)とし、音ごとの強度比を wとすれば、モデルの各正規分布の重み w^kは r^k -wと表される。これを用いて多重音基本周波数検出を行う。発音数推定プロセスでは、 1/妙らぎを考慮して = 1 /n (固定)と置き、基本周波数及び周波数パワー比の検出プロセスでは、 r^kの事前分布 p (r^k )を平均 l/n、分散 Vの正規分布と置くことで、最大事後確率推定を用いて周波数成分パワー比の推定を行う。 EMアルゴリズムにおける wと r^kの更新値はそれ k n

ぞれ式（20)、式（24)となる。

[数 21]

周波数成分パワー比パラメータを導入することで単一音と同等なスペクトル構造を有する和音を複数音として推定できる可能性がある。

産業上の利用可能性

[0049] 本発明は、音楽の音響信号形式 (CD, MD, .wavファイルなど)から、音楽演奏の形式 (MIDI信号， .midファイルなど)への変換 '編集'楽譜変換などのツールとして利用することができる。具体例を挙げると、着メロ作成、カラオケ作成、楽譜作成支援等に用いることができる。さら〖こ、本発明は、音源分離や音声認識へ応用することも可能である。

図面の簡単な説明

[0050] [図 1]多重音のスペクトルが、複数の調波構造の重ね合わせ力構成されることをと示す説明図である。

[図 2]観測された拡がりをもつスぺ外ルが正規分布で近似できることを説明する図である。

[図 3]EMアルゴリズムによる多重調波構造モデルの最尤推定を説明する図である。

[図 4]基本周波数 F力 371Hz, 441Hz, 556Hzの 3音による多重音スペクトルを示

0

す図であり、横軸は対数周波数、縦軸は振幅を表している。

[図 5]調波構造モデル数および基本周波数推定値の更新を示す図である。

[図 6]上図は観測されたスペクトルであって、横軸は対数周波数、縦軸は振幅であり、下図はが真の基本周波数 Fの 1/2であった場合の調波構造モデルを示す図であつて、横軸は対数周波数、縦軸は確率密度である。

圆 7]実験 1に用いた対象曲を示す表である。

[図 8]Kanonの楽譜の一部と同部分の基本周波数抽出結果を示す図である。

[図 9]上図は、話者 1人による音声信号における基本周波数 F検出結果を示す図で

0

あり、下図は、上図に対する参照基本周波数 Fパターンを示す図である。

0

[図 10]話者 1人による音声信号における実験結果を Cepstrum法の結果と併せて示す表である。

[図 11]上図は、話者 2人 (女性話者 2人）による音声信号における基本周波数 F検出

0 結果を示す図であり、下図は、上図に対する参照基本周波数 Fパターンを示す図で

0

ある。

[図 12]上図は、話者 2人 (男性話者 2人）による音声信号における基本周波数 F検出

0

ある。

[図 13]Cepstrum法の推定正解率を示す表である。

圆 14]本発明の手法の推定正解率を示す表である。

Claims

請求の範囲

[1] 多重音信号の観測スペクトルを、複数の調波構造モデルからなる多重調波構造モデルで推定することで多重音信号を解析する方法であって、

各調波構造モデルは、基本周波数成分に対応する一つの単峰形分布の代表値である基本周波数推定値と、該基本周波数推定値によって決定される他の単峰形分布の代表値とを有すると共に、調波構造モデルの各単峰形分布の代表値、重み、分散を多重調波構造モデルのモデルパラメータとし、

該方法は、モデルパラメータを推定するにあたり、所与の一つ又は複数の調波構造モデルに基づ、てモデル数を減少あるいは Zおよび増カロさせながら情報量基準を用いて適切な複数の調波構造モデルを選択し、選択された調波構造モデルの数を推定音源数とすることを特徴とする多重音信号解析方法。

[2] 請求項 1において、推定される音源数より多い数の調波構造モデルを用意し、尤度への寄与度が少ないと判断される調波構造モデルを削除しながら情報量基準を用いて適切な複数の調波構造モデルを選択することを特徴とする多重信号解析方法。

[3] 請求項 2にお、て、尤度への寄与度が少な!/、と判断される調波構造モデルは、調波構造モデルの調波構造ごとの重みが小さいモデルであることを特徴とする多重音信号解析方法。

[4] 請求項 3において、すべての調波構造モデルの中で重みが最小のモデル、あるいは、隣接する基本周波数推定値が所定の閾値より近接した場合には、重みが小さい方のモデル、のヽずれかに該当するモデルを削減することを特徴とする多重音信号解析方法。

[5] 請求項 1乃至 4いずれかの方法はさらに、選択された各調波構造モデルの基本周波数推定値から基本周波数を検出することを含むことを特徴とする多重音信号解析方法。

[6] 請求項 1乃至 4いずれかの方法はさらに、選択された各調波構造モデルの各単峰形分布の重みカゝら周波数成分パワー比を検出することを含むことを特徴とする多重音信号解析方法。

[7] 請求項 5, 6いずれかに記載の方法において、基本周波数あるいは zおよび周波数成分パワー比の検出は、

推定音源数の基となった各調波構造モデルにぉヽて、基本周波数推定値を構成する代表値を少なくとも整数倍および Zあるいは整数分の 1倍の値に順次置き換えながら少なくとも調波構造モデルの各単峰形分布の重みを推定するにあたり、情報量基準を用いて適切な置き換え値及び重みを選択し、選択された置き換え値ある、は Zおよび重みを、基本周波数あるいは Zおよび周波数成分パワー比の推定値とすることを特徴とする多重音信号解析方法。

[8] 多重音信号の観測スペクトルを、複数の調波構造モデルからなる多重調波構造モデルで推定することで多重音信号を解析する方法であって、

各調波構造モデルは、基本周波数成分に対応する一つの単峰形分布の代表値である基本周波数推定値と、該基本周波数推定値によって決定される他の単峰形分布の代表値とを有すると共に、調波構造モデルの各単峰形分布の代表値、重み、分散、調波モデル数を多重調波構造モデルのモデルパラメータとし、

該方法は、モデルパラメータのうち少なくとも基本周波数推定値を構成する代表値を推定するステップと、

推定された代表値を少なくとも整数倍および Zあるいは整数分の 1倍の値に順次置き換えながら少なくとも調波構造モデルの各単峰形分布の重みを推定するにあたり、情報量基準を用いて適切な置き換え値及び重みを選択するステップと、

を有し、選択された置き換え値あるいは Zおよび重みを基本周波数あるいは Zおよび周波数成分パワー比の推定値とすることを特徴とする多重音信号解析方法。

[9] 請求項 8にお、て、該方法は、多重音信号の音源数を推定するステップを含み、推定音源数の基となった各調波構造モデルにぉ、て基本周波数あるいは Zおよび周波数成分パワー比を検出するステップを適用することを特徴とする多重信号解析方法。

[10] 請求項 8において、多重音信号の音源数が既知であることを特徴とする多重音信号解析方法。

[11] 請求項 7乃至 10 、ずれかに記載の方法にぉ、て、

基本周波数推定値の置き換え値は、基本周波数推定値の少なくとも 1倍、 2倍、 3 倍、 2分の 1倍、 3分の 1倍、 3分の 2倍、 2分の 3倍からなる群から選択された 2つ以上の値であることを特徴とする多重音信号解析方法。

[12] 請求項 1乃至 11いずれか〖こおいて、単峰形分布は正規分布であることを特徴とする多重音信号解析方法。

[13] 請求項 1乃至 12いずれかにおいて、分布の代表値は平均であることを特徴とする多重音信号解析方法。

[14] 請求項 1乃至 13いずれかにおいて、分散を既知パラメータとして他のパラメータを推定することを特徴とする多重音信号解析方法。

[15] 請求項 1乃至 14いずれかにおいて、最尤推定法によってモデルパラメータの推定を行うことを特徴とする多重音信号解析法。

[16] 請求項 1乃至 14いずれかにおいて、各モデルパラメータに事前確率分布を想定して最大事後確率推定法 (MAP)によってモデルパラメータの推定を行うことを特徴とする多重信号解析法。

[17] 請求項 1乃至 14いずれかにおいて、変分ベイズ法によってモデルパラメータの推定を行うことを特徴とする多重信号解析法。

[18] 請求項 15において、前記情報量基準は AICであることを特徴とする多重信号解析法。

[19] 請求項 16において、前記情報量基準は BIC、 ABICあるいは MDLであることを特徴とする多重信号解析法。

[20] 請求項 1乃至 19いずれかにおいて、モデルパラメータの推定アルゴリズムは、 EMァルゴリズムであることを特徴とする多重音信号解析方法。

[21] 請求項 1乃至 20いずれかにおいて、モデルパラメータの推定は、観測スペクトルと多重調波構造モデル間の K-L情報量を評価尺度として行うことを特徴とする多重音信号解析方法。

[22] 請求項 1乃至 21いずれかにおいて、情報量基準を求める際の対数尤度は、観測スベクトルと多重調波構造モデル間の K-L情報量を最小化するパラメータを用いて算出されることを特徴とする多重音信号解析方法。

[23] 請求項 1乃至 22いずれかに記載の方法をコンピュータに実行させるためのコンビュタプログラム。

求項 1乃至 22いずれかに記載の方法をコンピュータに実行させるためのコンビュタプログラムを記録させた記録媒体。