(57)【要約】【課題】ノート指令入力信号パラメータと、ノート制御パラメータと、表現指令入力信号パラメータと、表現制御パラメータと、微細変動成分パラメータとを推定することができる。【解決手段】基本周波数抽出部１８は、音響信号の基本周波数軌跡を抽出する。初期状態系列推定部２２は、ノート指令入力信号及び表現指令入力信号の初期状態系列を推定する。信号分離部２４は、ノート成分、表現成分、及び微細変動成分の組からなる完全データの条件付き期待値及び完全データの自己相関の条件付き期待値を計算する。モデルパラメータ更新部２６は、目的関数をモデルパラメータに関して最大化するようにモデルパラメータを更新する。パラメータ収束判定部３６は、予め定められた終了条件を満たすまで、信号分離部２４及びモデルパラメータ更新部２６による処理を繰り返し、終了条件を満たすと判定された場合に、出力部１６が、モデルパラメータを出力する。【選択図】図３
2014026055.tif 000002
【特許請求の範囲】
【請求項１】
歌声を示す音響信号の基本周波数軌跡を、ノート指令入力信号とノート制御フィルタ特性を表すインパルス応答信号との畳みこみによって得られるノート成分と、表現指令入力信号と表現制御フィルタ特性を表すインパルス応答信号との畳みこみによって得られる表現成分と、微細変動成分との和で表したときの、前記ノート指令入力信号を表すモデルを構成するためのノート指令入力信号パラメータと、前記ノート制御フィルタ特性を表すモデルを構成するためのノート制御パラメータと、前記表現指令入力信号を表すモデルを構成するための表現指令入力信号パラメータと、前記表現制御フィルタ特性を表すモデルを構成するための表現制御パラメータと、前記微細変動成分を表すモデルを構成するための微細変動成分パラメータとを推定する音響信号分析方法であって、入力された歌唱旋律の各音符の音高を示す音高目標値時系列、及び前記ノート指令入力信号パラメータと前記ノート制御パラメータと前記表現指令入力信号パラメータと前記表現制御パラメータと前記微細変動成分パラメータとからなるモデルパラメータの初期値を記憶したパラメータ初期値データベースと、基本周波数抽出手段と、初期状態系列推定手段と、信号分離手段と、モデルパラメータ更新手段と、収束判定手段と、パラメータ出力手段とを含む音響信号分析装置における音響信号分析方法であって、
前記基本周波数抽出手段によって、歌唱者の歌声を示す音響信号を入力として、前記音響信号の基本周波数軌跡を抽出し、
前記初期状態系列推定手段によって、前記音高目標値時系列と前記基本周波数軌跡とに基づいて、前記ノート指令入力信号及び前記表現指令入力信号の各時刻の状態を示す初期状態系列を推定し、
前記信号分離手段によって、前記音高目標値時系列、前記基本周波数軌跡、及び前記モデルパラメータの初期値または前回更新された前記モデルパラメータに基づいて、前記基本周波数軌跡と、前記モデルパラメータとが与えられたときの、前記ノート成分、前記表現成分、及び前記微細変動成分の組からなる完全データの条件付き期待値及び前記完全データの自己相関の条件付き期待値を計算し、
前記モデルパラメータ更新手段によって、前記信号分離手段によって計算された前記完全データの条件付き期待値及び前記完全データの自己相関の条件付き期待値、並びに前記モデルパラメータの事前確率を用いて表わされる、前記完全データの対数尤度関数の期待値を表わす関数を、目的関数とし、前記音高目標値時系列、前記基本周波数軌跡、及び前記モデルパラメータの初期値または前回更新された前記モデルパラメータに基づいて、前記目的関数を前記モデルパラメータに関して最大化するように前記モデルパラメータを更新し、
前記収束判定手段によって、予め定められた終了条件を満たすまで、前記信号分離手段による計算及び前記更新手段による更新を繰り返し、
前記パラメータ出力手段によって、前記パラメータ収束判定手段により前記終了条件を満たすと判定された場合に、前記モデルパラメータを出力する
音響信号分析方法。
【請求項２】
前記ノート指令入力信号は階段状の信号であり、各時刻の状態を示す状態系列パラメータと、各音符の音高からの推移を示す音高シフトパラメータとからなる前記ノート指令入力信号パラメータを用いてモデル化される請求項１記載の音響信号分析方法。
【請求項３】
前記表現指令入力信号は矩形状のパルス信号であり、各時刻の状態を示す状態系列パラメータと、各音符の各状態における表現指令の大きさを表わす表現振幅パラメータとからなる前記表現指令入力信号パラメータを用いて、隠れマルコフモデルによって確率モデル化される請求項１又は２記載の音響信号分析方法。
【請求項４】
前記微細変動成分は、前記微細変動成分パラメータからなるガウス性白色雑音として確率モデル化される請求項１〜請求項３の何れか１項記載の音響信号分析方法。
【請求項５】
前記信号生成系のノート制御フィルタ特性は、線形２次系システムによって導出されるフィルタで表され、前記ノート制御パラメータは、線形２次系システムの固有周波数に反比例するパラメータと、線形２次系システムの減衰率パラメータである請求項１〜請求項４の何れか１項記載の音響信号分析方法。
【請求項６】
前記信号生成系の表現制御フィルタ特性は、臨界制動２次系システムによって導出されるフィルタで表され、前記表現制御パラメータは、臨界制動２次系システムの固有周波数に反比例するパラメータである請求項１〜請求項５の何れか１項記載の音響信号分析方法。
【請求項７】
歌声を示す音響信号の基本周波数軌跡を、ノート指令入力信号とノート制御フィルタ特性を表すインパルス応答信号との畳みこみによって得られるノート成分と、表現指令入力信号と表現制御フィルタ特性を表すインパルス応答信号との畳みこみによって得られる表現成分と、微細変動成分との和で表したときの、前記ノート指令入力信号を表すモデルを構成するためのノート指令入力信号パラメータと、前記ノート制御フィルタ特性を表すモデルを構成するためのノート制御パラメータと、前記表現指令入力信号を表すモデルを構成するための表現指令入力信号パラメータと、前記表現制御フィルタ特性を表すモデルを構成するための表現制御パラメータと、前記微細変動成分を表すモデルを構成するための微細変動成分パラメータとを推定する音響信号分析装置であって、
歌唱者の歌声を示す音響信号を入力として、前記音響信号の基本周波数軌跡を抽出する基本周波数抽出手段と、
入力された歌唱旋律の各音符の音高を示す音高目標値時系列、及び前記ノート指令入力信号パラメータと前記ノート制御パラメータと前記表現指令入力信号パラメータと前記表現制御パラメータと前記微細変動成分パラメータとからなるモデルパラメータの初期値を記憶したパラメータ初期値データベースと、
前記音高目標値時系列と前記基本周波数軌跡とに基づいて、前記ノート指令入力信号及び前記表現指令入力信号の各時刻の状態を示す初期状態系列を推定する初期状態系列推定手段と、
前記音高目標値時系列、前記基本周波数軌跡、及び前記モデルパラメータの初期値または前回更新された前記モデルパラメータに基づいて、前記基本周波数軌跡と、前記モデルパラメータとが与えられたときの、前記ノート成分、前記表現成分、及び前記微細変動成分の組からなる完全データの条件付き期待値及び前記完全データの自己相関の条件付き期待値を計算する信号分離手段と、
前記信号分離手段によって計算された前記完全データの条件付き期待値及び前記完全データの自己相関の条件付き期待値、並びに前記モデルパラメータの事前確率を用いて表わされる、前記完全データの対数尤度関数の期待値を表わす関数を、目的関数とし、前記音高目標値時系列、前記基本周波数軌跡、及び前記モデルパラメータの初期値または前回更新された前記モデルパラメータに基づいて、前記目的関数を前記モデルパラメータに関して最大化するように前記モデルパラメータを更新するモデルパラメータ更新手段と、
予め定められた終了条件を満たすまで、前記信号分離手段による計算及び前記更新手段による更新を繰り返し行う収束判定手段と、
前記パラメータ収束判定部により前記終了条件を満たすと判定された場合に、前記モデルパラメータを出力するパラメータ出力手段と、
を含む音響信号分析装置。
【請求項８】
コンピュータに、請求項１〜請求項６の何れか１項記載の音響信号分析方法の各ステップを実行させるためのプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、音響信号分析方法、装置、及びプログラムに係り、歌声を示す音響信号の基本周波数軌跡を表わすパラメータを推定する音響信号分析方法、装置、及びプログラムに関する。
【背景技術】
【０００２】
楽曲のメロディを歌った歌声のＦ０軌跡には、図６に示すように、そのメロディを構成する楽譜の音符の並びだけでなく、楽譜に記載されない、“楽譜から逸脱した”動的変動成分が含まれる。これらは、発声器官の物理的制約に起因する成分（特に、オーバーシュートやプレパレーション、微細変動成分）と意図的表現による成分（特に、ビブラートやポルタメント）からなると考えられ、知覚的には、前者は人間らしさ・自然性に関係し、後者は巧拙感に関係することがわかってきている。さらに、後者は意図して意図通りにずれた場合と、意図通りに制御できなかったずれとに分かれ、習熟度に関連すると考えられる。このような物理的制約もしくは意図的表現による楽譜から逸脱した動的変動成分をＦ０軌跡から楽譜と分離して抽出できれば、歌唱者の歌い方や個性、癖を学習することにつながり、現在盛んに研究される歌声合成や歌声変換などへの応用が期待できる。例えば、ある歌声を別の歌唱者の歌い方に変換して合成することが可能となるだろう。従来と違って、事前に歌唱者の歌い方が学習されるため、どんなメロディにもその歌い方を転写できることを特徴とする。
【０００３】
従来、線形２次系システムを利用して歌声の基本周波数に含まれる動的変動成分を制御するモデルが提案されている（非特許文献１、非特許文献２、非特許文献３）。これらの研究では、日本語の話声の基本周波数パターンを表現する藤崎モデルが参考にされている。藤崎モデルは、臨界制動２次系のインパルス応答とステップ応答を利用して、日本語の句頭から句末に向けて緩やかに下降するフレーズ成分と、語句に対応して急激に上昇下降するアクセント成分を表現し、これらを重畳することで、基本周波数軌跡を記述する。ただし、歌声の旋律に伴った急激な基本周波数の上昇・下降の制御及び、ビブラートのような周期的な振動は、臨界制動系では表現できない。そのため、歌声のＦ０制御モデルでは２次系の伝達関数
【０００４】
【数１】
2014026055.tif 000003
【０００５】
における減衰率ζを調整することによって、指数減衰（ζ＞１）、減衰振動（０＜ζ＜１、オーバーシュートに対応する）、臨界制動（ζ＝１）、定常振動（ζ＝０、ビブラートに対応する）からなる様々な振動現象を表現する。また、楽譜の音符列を表す階段状信号に上記（１）式のインパルス応答を部分的に畳み込んで得られる基本周波数軌跡を利用して、表情豊かな歌声合成音を実現している（非特許文献３）。しかしながら、これらの従来技術では、制御パラメータ（減衰率ζと固有周波数Ω）が手作業あるいは規則に基づいて決定されるものであった。
【０００６】
これに対し、観測されるＦ０軌跡からモデルパラメータを推定する逆問題の解法が提案されている（非特許文献４）。ここでは、Ｆ０軌跡の生成過程を図７の下部のように定義し、ノートの変化時点を始点終点と考え、それによって分割されるセグメントごとに、信号生成系の各パラメータを反復的に推定する信号解析方法を提供している。また、同じ枠組みで藤崎モデルのパラメータ推定法も提供されている（非特許文献５）。
【先行技術文献】
【非特許文献】
【０００７】
【非特許文献１】N. Minematsu, B. Matsuoka, and K. Hirose, “ Prosodic Modeling of Nagauta Singing and Its Evaluation, ” in Proc. SpeechProsody 2004, pp. 487−490, Mar. 2004.
【非特許文献２】T. Saitou, M. Unoki, and M. Akagi,“ Development of an F0 control Model Based on F0 Dynamic Characteristics for Singing-Voice Synthesis, ” Speech Communication, vol.46, pp. 405−417, 2005.
【非特許文献３】T. Saitou, M. Goto, M. Unoki, and M. Akagi, “ Speech-To-Singing Synthesis: Converting Speaking Voices to Singing Voices by Controlling Acoustic Features Unique to Singing Voices, ”in Proc. WASPAA 2007, pp. 215−218, Oct. 2007.
【非特許文献４】Y. Ohishi, H. Kameoka, D. Mochihashi, H. Nagano, K. Kashino, “ Statistical Modeling of F0 Dynamics in Singing Voices Based on Gaussian Processes with Multiple Oscillation Bases, ”in Proc. INTERSPEECH 2010, pp. 2598−2601, Sept. 2010.
【非特許文献５】H. Kameoka, J. Le Roux, Y. Ohishi,“A statistical model of speech F0 contours, ”ISCA Tutorial and Research Workshop on Statistical And Perceptual Audition (SAPA 2010), pp. 43-48, Sept. 2010.
【発明の概要】
【発明が解決しようとする課題】
【０００８】
しかしながら、上記の非特許文献４に記載の技術では、モデルパラメータを推定するために、Ｆ０軌跡をあらかじめ、音符区間に対応するセグメントに分割する必要があった（手動による分割、もしくは楽譜とのビタビアライメントによる分割）。また、ビブラートや微細変動成分のすべてがガウス性白色雑音としてモデル化されていたため、歌唱者の表現意図による動的変動成分（ビブラートなどの表現成分）を微細変動成分と分離して、明確に特徴付けられなかった。
【０００９】
本発明では、上記の事情を考慮してなされたものであり、歌声を示す音響信号の基本周波数軌跡を、ノート成分と、表現成分と、微細変動成分との和で表したときの、ノート指令入力信号パラメータと、ノート制御パラメータと、表現指令入力信号パラメータと、表現制御パラメータと、微細変動成分パラメータとを推定することができる音響信号分析方法、装置及びプログラムを提供することを目的とする。
【課題を解決するための手段】
【００１０】
上記目的を達成するために、第１の発明の音響信号分析方法は、歌声を示す音響信号の基本周波数軌跡を、ノート指令入力信号とノート制御フィルタ特性を表すインパルス応答信号との畳みこみによって得られるノート成分と、表現指令入力信号と表現制御フィルタ特性を表すインパルス応答信号との畳みこみによって得られる表現成分と、微細変動成分との和で表したときの、前記ノート指令入力信号を表すモデルを構成するためのノート指令入力信号パラメータと、前記ノート制御フィルタ特性を表すモデルを構成するためのノート制御パラメータと、前記表現指令入力信号を表すモデルを構成するための表現指令入力信号パラメータと、前記表現制御フィルタ特性を表すモデルを構成するための表現制御パラメータと、前記微細変動成分を表すモデルを構成するための微細変動成分パラメータとを推定する音響信号分析方法であって、入力された歌唱旋律の各音符の音高を示す音高目標値時系列、及び前記ノート指令入力信号パラメータと前記ノート制御パラメータと前記表現指令入力信号パラメータと前記表現制御パラメータと前記微細変動成分パラメータとからなるモデルパラメータの初期値を記憶したパラメータ初期値データベースと、基本周波数抽出手段と、初期状態系列推定手段と、信号分離手段と、モデルパラメータ更新手段と、収束判定手段と、パラメータ出力手段とを含む音響信号分析装置における音響信号分析方法であって、前記基本周波数抽出手段によって、歌唱者の歌声を示す音響信号を入力として、前記音響信号の基本周波数軌跡を抽出し、前記初期状態系列推定手段によって、前記音高目標値時系列と前記基本周波数軌跡とに基づいて、前記ノート指令入力信号及び前記表現指令入力信号の各時刻の状態を示す初期状態系列を推定し、前記信号分離手段によって、前記音高目標値時系列、前記基本周波数軌跡、及び前記モデルパラメータの初期値または前回更新された前記モデルパラメータに基づいて、前記基本周波数軌跡と、前記モデルパラメータとが与えられたときの、前記ノート成分、前記表現成分、及び前記微細変動成分の組からなる完全データの条件付き期待値及び前記完全データの自己相関の条件付き期待値を計算し、前記モデルパラメータ更新手段によって、前記信号分離手段によって計算された前記完全データの条件付き期待値及び前記完全データの自己相関の条件付き期待値、並びに前記モデルパラメータの事前確率を用いて表わされる、前記完全データの対数尤度関数の期待値を表わす関数を、目的関数とし、前記音高目標値時系列、前記基本周波数軌跡、及び前記モデルパラメータの初期値または前回更新された前記モデルパラメータに基づいて、前記目的関数を前記モデルパラメータに関して最大化するように前記モデルパラメータを更新し、前記収束判定手段によって、予め定められた終了条件を満たすまで、前記信号分離手段による計算及び前記更新手段による更新を繰り返し、前記パラメータ出力手段によって、前記パラメータ収束判定手段により前記終了条件を満たすと判定された場合に、前記モデルパラメータを出力する。
【００１１】
第２の発明の音響信号分析装置は、歌声を示す音響信号の基本周波数軌跡を、ノート指令入力信号とノート制御フィルタ特性を表すインパルス応答信号との畳みこみによって得られるノート成分と、表現指令入力信号と表現制御フィルタ特性を表すインパルス応答信号との畳みこみによって得られる表現成分と、微細変動成分との和で表したときの、前記ノート指令入力信号を表すモデルを構成するためのノート指令入力信号パラメータと、前記ノート制御フィルタ特性を表すモデルを構成するためのノート制御パラメータと、前記表現指令入力信号を表すモデルを構成するための表現指令入力信号パラメータと、前記表現制御フィルタ特性を表すモデルを構成するための表現制御パラメータと、前記微細変動成分を表すモデルを構成するための微細変動成分パラメータとを推定する音響信号分析装置であって、歌唱者の歌声を示す音響信号を入力として、前記音響信号の基本周波数軌跡を抽出する基本周波数抽出手段と、入力された歌唱旋律の各音符の音高を示す音高目標値時系列、及び前記ノート指令入力信号パラメータと前記ノート制御パラメータと前記表現指令入力信号パラメータと前記表現制御パラメータと前記微細変動成分パラメータとからなるモデルパラメータの初期値を記憶したパラメータ初期値データベースと、前記音高目標値時系列と前記基本周波数軌跡とに基づいて、前記ノート指令入力信号及び前記表現指令入力信号の各時刻の状態を示す初期状態系列を推定する初期状態系列推定手段と、前記音高目標値時系列、前記基本周波数軌跡、及び前記モデルパラメータの初期値または前回更新された前記モデルパラメータに基づいて、前記基本周波数軌跡と、前記モデルパラメータとが与えられたときの、前記ノート成分、前記表現成分、及び前記微細変動成分の組からなる完全データの条件付き期待値及び前記完全データの自己相関の条件付き期待値を計算する信号分離手段と、前記信号分離手段によって計算された前記完全データの条件付き期待値及び前記完全データの自己相関の条件付き期待値、並びに前記モデルパラメータの事前確率を用いて表わされる、前記完全データの対数尤度関数の期待値を表わす関数を、目的関数とし、前記音高目標値時系列、前記基本周波数軌跡、及び前記モデルパラメータの初期値または前回更新された前記モデルパラメータに基づいて、前記目的関数を前記モデルパラメータに関して最大化するように前記モデルパラメータを更新するモデルパラメータ更新手段と、予め定められた終了条件を満たすまで、前記信号分離手段による計算及び前記更新手段による更新を繰り返し行う収束判定手段と、前記パラメータ収束判定部により前記終了条件を満たすと判定された場合に、前記モデルパラメータを出力するパラメータ出力手段と、を含んで構成されている。
【００１２】
第１の発明及び第２の発明によれば、前記基本周波数抽出手段によって、歌唱者の歌声を示す音響信号を入力として、前記音響信号の基本周波数軌跡を抽出し、前記初期状態系列推定手段によって、前記音高目標値時系列と前記基本周波数軌跡とに基づいて、前記ノート指令入力信号及び前記表現指令入力信号の各時刻の状態を示す初期状態系列を推定し、前記信号分離手段によって、前記音高目標値時系列、前記基本周波数軌跡、及び前記モデルパラメータの初期値または前回更新された前記モデルパラメータに基づいて、前記基本周波数軌跡と、前記モデルパラメータとが与えられたときの、前記ノート成分、前記表現成分、及び前記微細変動成分の組からなる完全データの条件付き期待値及び前記完全データの自己相関の条件付き期待値を計算し、前記モデルパラメータ更新手段によって、前記信号分離手段によって計算された前記完全データの条件付き期待値及び前記完全データの自己相関の条件付き期待値、並びに前記モデルパラメータの事前確率を用いて表わされる、前記完全データの対数尤度関数の期待値を表わす関数を、目的関数とし、前記音高目標値時系列、前記基本周波数軌跡、及び前記モデルパラメータの初期値または前回更新された前記モデルパラメータに基づいて、前記目的関数を前記モデルパラメータに関して最大化するように前記モデルパラメータを更新し、前記収束判定手段によって、予め定められた終了条件を満たすまで、前記信号分離手段による計算及び前記更新手段による更新を繰り返し、前記パラメータ出力手段によって、前記パラメータ収束判定手段により前記終了条件を満たすと判定された場合に、前記モデルパラメータを出力する。
【００１３】
このように、歌声を示す音響信号の基本周波数軌跡を、ノート成分と、表現成分と、微細変動成分との和で表したときの、ノート指令入力信号パラメータと、ノート制御パラメータと、表現指令入力信号パラメータと、表現制御パラメータと、微細変動成分パラメータとを推定することができる。
【００１４】
また、ノート指令入力信号は階段状の信号であり、各時刻の状態を示す状態系列パラメータと、各音符の音高からの推移を示す音高シフトパラメータとからなるノート指令入力信号パラメータを用いてモデル化される。
【００１５】
また、表現指令入力信号は矩形状のパルス信号であり、各時刻の状態を示す状態系列パラメータと、各音符の各状態における表現指令の大きさを表わす表現振幅パラメータとからなる表現指令入力信号パラメータを用いて、隠れマルコフモデルによって確率モデル化される。
【００１６】
また、微細変動成分は、微細変動成分パラメータからなるガウス性白色雑音として確率モデル化される。
【００１７】
また、信号生成系のノート制御フィルタ特性は、線形２次系システムによって導出されるフィルタで表され、前記ノート制御パラメータは、線形２次系システムの固有周波数に反比例するパラメータと、線形２次系システムの減衰率パラメータである。
【００１８】
また、信号生成系の表現制御フィルタ特性は、臨界制動２次系システムによって導出されるフィルタで表され、前記表現制御パラメータは、臨界制動２次系システムの固有周波数に反比例するパラメータである。
【００１９】
本発明のプログラムは、コンピュータに音響信号分析方法の各ステップを実行させるためのプログラムである。
【発明の効果】
【００２０】
以上説明したように、本発明の音響信号分析方法、装置、及びプログラムによれば、歌声を示す音響信号の基本周波数軌跡を、ノート成分と、表現成分と、微細変動成分との和で表したときの、ノート指令入力信号パラメータと、ノート制御パラメータと、表現指令入力信号パラメータと、表現制御パラメータと、微細変動成分パラメータとを推定することができる。
【図面の簡単な説明】
【００２１】
【図１】提案モデルの概略図である。
【図２】隠れマルコフモデルによる指令入力信号の統計モデル化を示す図である。
【図３】本実施の形態の音響信号分析装置の機能的構成を示すブロック図である。
【図４】本実施の形態の音響信号分析装置におけるパラメータ推定処理ルーチンの内容を示すフローチャートである。
【図５】テストデータを示す図である。
【図６】歌声の基本周波数軌跡（Ｆ０）を示す図である。
【図７】従来技術におけるＦ０生成過程を示す図である。
【発明を実施するための形態】
【００２２】
以下、図面を参照して本発明の実施の形態を詳細に説明する。本実施の形態では、後述するパラメータ推定アルゴリズムを用いて信号解析を行う構成である。
【００２３】
＜発明の概要＞
本発明では、歌唱する楽譜の音符列を表すノート指令入力信号と表現指令入力信号によって駆動される歌声Ｆ０軌跡の生成過程モデルを提案し、歌声の基本周波数（Ｆ０）軌跡から、オーバーシュートやポルタメントを表現するノート制御機構（フィルタ）のパラメータ（ノート制御パラメータと呼ぶ）およびノート指令入力信号を構成するパラメータ（ノート指令入力信号パラメータと呼ぶ）、ビブラートなどの表現意図を表す表現制御機構（フィルタ）のパラメータ（表現制御パラメータと呼ぶ）および表現指令入力信号を構成するパラメータ（表現指令入力信号パラメータと呼ぶ）、そして微細変動成分を表現するパラメータ（微細変動成分パラメータと呼ぶ）を抽出する。ここで、ノート指令入力信号は、図１に示す通り、楽譜に記載される音符の並び（階段状に音階が変化する）を表現する。一方で、表現指令入力信号は歌唱者の音楽的な表現意図（ノート指令入力信号に比べて細かいステップ状の指令）を表現する。このモデルは、対数スケールのＦ０軌跡ｙ（ｔ）（ここで、ｔは時間を表す）が３つの成分の重ね合わせで表現されることを想定する。図１に示されるノート成分と表現成分は先に述べた指令信号によって駆動されるノート制御機構、表現制御機構の出力である。ノート制御機構と表現制御機構はある種のフィルタに相当し、線形２次系によって表現される。ノート制御機構は、オーバーシュート、ポルタメントのようなノート（音符）の立ち上がり方を線形２次系によって制御する。表現制御機構は、ビブラートなどの音楽的な意図表現を線形２次系によって制御する。ただし、表現制御機構は、線形２次系に含まれる１つの系、「臨界制動２次系（上記（１）式のζ＝１の場合）」で構成する（非特許文献３）。最後に、図１に示される微細変動成分は、１０Ｈｚ以上の不規則な振動成分を想定する（非特許文献３）。本発明では、このようなＦ０軌跡の生成過程を離散時間領域で表現し、統計的手法に基づいてＦ０軌跡生成系のモデルパラメータを推定するためにフレームワークを提供する。
＜原理＞
次に、Ｆ０軌跡生成系のモデルパラメータを推定する原理について説明する。
＜提案モデルの離散時間表現＞
まず、本実施形態の提案モデルの離散時間表現について説明する。
【００２４】
連続時間領域で表現されるノート制御機構、表現制御機構の２次系の伝達関数の離散時間表現を得るために、後退差分変換を利用する。後退差分変換は、時間微分演算子ｓをｚ領域における後退差分演算子ｓ≒（１−ｚ-1）／ｔ0に置き換える変換であり、（ｔ0は離散時間表現におけるサンプリング周期とする）、この変換によりノート制御機構の逆システムの伝達関数Ｈn-1（ｓ）はｚ領域で、
【００２５】
【数２】
2014026055.tif 000004
と書くことができる。ただし、
【００２６】
【数３】
2014026055.tif 000005
および、φ、ψは、φ＝１／（Ωｔ0）、ψ＝ζと表現される。ここで、ｋを離散時刻インデックスとし、ノート指令入力信号およびノート成分の離散時間表現をそれぞれｕn［ｋ］、ｙｎ［ｋ］とすると、ｙｎ［ｋ］は、ノート制御パラメータφ、ψによって特性が決まる拘束つき全極モデルからの出力
【００２７】
【数４】
2014026055.tif 000006
とみなすことができる。同様に、表現指令入力信号ｕe［ｋ］と表現成分ｙe［ｋ］の関係も
【００２８】
【数５】
2014026055.tif 000007
と書くことができる。ただし、ｂ2＝ξ2,ｂ1＝−２ξ（１＋ξ）,ｂ0＝１＋２ξ＋ξ2であり、表現制御パラメータξはξ＝１／（Ωｔ0）と表現される。微細変動成分ｙf（ｔ）の離散時間表現をｙf［ｋ］とすると、提案モデルによる歌声Ｆ０軌跡の離散時間表現は、これら３つの成分の和
【００２９】
【数６】
2014026055.tif 000008
で与えられる。
【００３０】
＜歌声のＦ０生成モデルの統計モデル化＞
次に、歌声のＦ０生成モデルの統計モデル化について説明する。
【００３１】
ノート指令入力信号と表現指令入力信号は図１に示す通り、それぞれ、楽譜に記載されるメロディの音符の並び、および歌唱者の音楽的な表現意図を表すものと想定する。この想定を指令入力信号に組み込むために、隠れマルコフモデル（ＨＭＭ）を利用して、ｕn［ｋ］とｕe［ｋ］を確率モデル化する。まず、ｏ［ｋ］：＝（ｕn［ｋ］、ｕe［ｋ］）Ｔを、
【００３２】
【数７】
2014026055.tif 000009
のように正規分布する確率変数と見なし、平均ν［ｋ］が、図２のようなノート指令入力信号と表現指令入力信号の組み合わせの状態に関する状態遷移に伴って変化するモデルを考える。これはＨＭＭに他ならず、このようにｏ［ｋ］をＨＭＭでモデル化したことにより、状態遷移の経路制限（状態遷移確率の設定）を通して、ν［ｋ］に制約を与えることが可能となる。
【００３３】
具体的には、図２に示す通り、このＨＭＭはＩ×Ｊ個の状態からなる。これらの各状態では、μn［ｋ］はＡn(i)＋ｄiの値をとる。ここで、Ａn(i)は楽譜から与えられる値であり，楽譜に記載されるメロディのi番目の音符の音高を表し、ｄiはその音高からの推移（音高シフトパラメータと呼ぶ）、Ｉは歌唱するメロディに含まれる音符の総数に相当する。したがって、μn［ｋ］はI個の音符からなる階段状信号を表現する。一方、ｉ番目の音符は、Ｓi,1を通らずして、状態Ｓi,jから別の状態Ｓi,j'（ｊ≠ｊ’，２≦ｊ≦Ｊ，２≦ｊ’≦Ｊ） へ直接に遷移できない。この制約によって、歌唱者の表現意図を表す表現指令入力信号μe［ｋ］は図２の下段に示すような矩形状のパルス信号を表現する。このとき、表現指令の大きさを｛Ｂe(i,j)｝i=1,j=1I,Jと表し、Ｊは表現指令の大きさ０を含む各音符における指令の種類数に相当する（図２では、Ｊ＝３とする）。以上より、指令入力信号を生成するＨＭＭの構成は以下のとおりである。
【００３４】
【数８】
2014026055.tif 000010
【００３５】
簡単のため状態遷移確率φn',nは定数とすると、指令入力モデルにおいて推定すべきパラメータは、ＨＭＭのどういった経路を辿るか表現する状態系列｛ｓk｝k=1K、音高シフトパラメータ｛ｄi｝i=1I、表現指令の大きさパラメータ｛Ｂe(i,j)｝i=1,j=1I,J、入力指令信号の分散パラメータσn2、σe2であり、これらをまとめて
【００３６】
【数９】
2014026055.tif 000011
と記述する。なお、音高シフトパラメータ｛ｄi｝i=1I、分散パラメータσn2が、ノート指令入力信号パラメータであり、表現指令の大きさパラメータ｛Ｂe(i,j)｝i=1,j=1I,J、分散パラメータσe2が、表現指令入力信号パラメータである。また、平均値系列｛μn［ｋ］｝k=1Kおよび｛μe［ｋ］｝k=1Kは、状態系列｛ｓk｝k=1Kが与えられたもとで、
【００３７】
【数１０】
2014026055.tif 000012
で与えられる。
【００３８】
指令入力信号モデルに基づいて、ｙ ＝ （ｙ［１］，．．．，ｙ［Ｋ］）Tの確率密度関数を導く。上記（７）式より、ｕｎ：＝ （ｕｎ［１］，．．．，ｕｎ［Ｋ］）T, ｕe：＝ （ｕe［１］，．．．，ｕe［Ｋ］）T, μｎ：＝ （μｎ［１］ ，．．．， μｎ［Ｋ］）T, μe：＝ （μe［１］ ，．．．， μe［Ｋ］）Tとすると
【００３９】
【数１１】
2014026055.tif 000013
が言える。ここで、Ｎ（μn,Σn）は、平均μn、分散Σnのガウス分布を表わす。ＩはＫ×Ｋの単位行列を表す。ノート成分ｙｎ：＝ （ｙｎ［１］，．．．，ｙｎ［Ｋ］）Tとノート指令入力信号ｕｎの関係、および表現成分ｙe：＝ （ｙe［１］，．．．，ｙe［Ｋ］）Tと表現指令入力信号ｕeの関係は、
【００４０】
【数１２】
2014026055.tif 000014
と置くと、それぞれ、
【００４１】
【数１３】
2014026055.tif 000015
と表現できる。ここでは、ノート制御パラメータφ、ψがノートごとに異なるものと想定し、｛φ(i)、ψ(i)｝i=1Iをパラメータとして拡張する（非特許文献３）。したがって、ノート成分と表現成分の確率密度関数は、
【００４２】
【数１４】
2014026055.tif 000016
が導出される。微細変動成分ｙf［Ｋ］はガウス性白色雑音を想定し、
【００４３】
【数１５】
2014026055.tif 000017
と表現する。なお、分散パラメータσf2が、微細変動成分パラメータである。
【００４４】
仮定より、ｙn, ｙe, ｙfは独立なので、モデルパラメータΘ：＝｛θu, ｛φ(i)、ψ(i)｝i=1I,ξ,σf2｝が与えられたもとでのＦ０軌跡ｙ＝ ｙn＋ｙe＋ｙfの確率密度関数は上記（１３）式、（１４）式と（１５）式より、
【００４５】
【数１６】
2014026055.tif 000018
で与えられる。以上より、
【００４６】
【数１７】
2014026055.tif 000019
が歌声Ｆ０軌跡ｙが与えられたときのモデルパラメータΘの尤度関数である。
【００４７】
Θの事前確率については、各パラメータは互いに独立で、パラメータσn2,σe2,σf2は一様に分布すると仮定し、
【００４８】
【数１８】
2014026055.tif 000020
とする。
【００４９】
＜パラメータ推定アルゴリズム＞
次に、パラメータ推定アルゴリズムについて説明する。
【００５０】
ｙが与えられた下で、Θの事後確率Ｐ（Θ｜ｙ)を最大化する問題は解析的に解くことはできないが、ｘ：＝（ｙnT,ｙeT,ｙfT)Tを完全データと見なすことで、ＥＭアルゴリズムによる不完全データ問題に帰着できる。不完全データと完全データの関係はｙ＝Ｈxとなり、ここでＨ：＝［Ｉ Ｉ Ｉ］とする。この場合、完全データの対数尤度は、
【００５１】
【数１９】
2014026055.tif 000021
で与えられる。このとき、完全データの対数尤度の期待値を表わすＱ関数Ｑ(Θ,Θ´)は、
【００５２】
【数２０】
2014026055.tif 000022
となる。ここで、ｔｒ（・）は行列のトレースを表し、Ｅ［ｘ｜ｙ；Θ´］とＥ［ｘｘT｜ｙ；Θ´］は、完全データの条件付き期待値、完全データの自己相関の条件付き期待値であり、条件付きガウス分布の性質より、
【００５３】
【数２１】
2014026055.tif 000023
と書ける。Ｅステップでは、直前のステップで更新されたパラメータをΘ´に代入し、上記に基づいてＥ［ｘ｜ｙ； Θ´］とＥ［ｘｘT｜ｙ； Θ´］が算出される。ｙn,ｙe,ｙfに対応するように、Ｅ［ｘ｜ｙ； Θ´］及びＥ［ｘｘT｜ｙ； Θ´］を
【００５４】
【数２２】
2014026055.tif 000024
のように区分表現すると、Ｑ関数は
【００５５】
【数２３】
2014026055.tif 000025
と書き直せて、これを用いて各パラメータについてＭステップの更新式を求めることができる。
【００５６】
＜状態系列＞
Ｑ関数の中でｓ：＝｛ｓk｝k=１Kに関する項は
【００５７】
【数２４】
2014026055.tif 000026
となる。ここで、ｏ［ｋ］：＝（［Ａ−ｘｎ］k,［Ｂ−ｘe]k)Tであり、［・］kはベクトルのｋ番目の要素を表す。これを最大化する状態系列｛ｓk｝k=1Kは動的計画法により効率的に解くことができる。まず、状態Ｓ1,1について
【００５８】
【数２５】
2014026055.tif 000027
とおくと、ｋ＝２，３，．．．，Ｋについて逐次的にδk（Ｓi,j）を
【００５９】
【数２６】
2014026055.tif 000028
により計算できる。各ステップで選択される状態番号
【００６０】
【数２７】
2014026055.tif 000029
を記憶しておくことで、ｋ＝Ｋまで到着後にｓk-1＝ψk（sk）（ｋ＝Ｋ，Ｋ−１，．．．，２）により選択された状態番号を辿っていくと最適経路ｓ1，．．．， ｓkを得る。
【００６１】
＜ノート制御パラメータ＞
φ(i)とψ(i)に関する事前分布をφ(i)〜Ｎ（μφ,σ2φ）、ψ(i)〜Ｎ（μψ,σ2ψ）とする。Ｑ関数の中でφ(i)とψ(i)に関係する項は、
【００６２】
【数２８】
2014026055.tif 000030
となる。ここで｜Ｔ｜は集合Ｔの要素数を表す。また、［μ］Ｔは、集合Ｔの要素をインデックスとして、μからそのインデックスに相当する要素を取り出した部分ベクトルを表す。今、
【００６３】
【数２９】
2014026055.tif 000031
として、上記(4)式から、Ａ(i)は、
【００６４】
【数３０】
2014026055.tif 000032
と表現される。ここで、［Ｕ］Ｔは集合Ｔの要素をインデックスとして、行列Ｕからそのインデックスに相当する行ベクトルを取り出して構成される部分行列を意味する。ニュートン・ラフソン法を利用して、Ｉ2(φ(i)、ψ(i))を最大化するφ(i)とψ(i)が数値的に導出される。
【００６５】
＜表現制御パラメータ＞
ξに関する事前分布をξ〜Ｎ（μξ,σξ2）とする。Ｑ関数の中でξに関係する項は、
【００６６】
【数３１】
2014026055.tif 000033
となる。ニュートン・ラフソン法を利用して、Ｉ3（ξ）を最大化するξが数値的に導出される。
【００６７】
＜その他のパラメータ＞
ｄiとＢe(i,j)に関して、それぞれ事前分布をｄi〜Ｎ（０,σｄ2）と
【００６８】
【数３２】
2014026055.tif 000034
とする。残されたパラメータの更新式は、
【００６９】
【数３３】
2014026055.tif 000035
と導出される。ＥＭアルゴリズムに基づくパラメータ推定アルゴリズムをまとめると、
【００７０】
初期化：パラメータΘの初期化
E-step：Ｅ［ｘ｜ｙ；Θ´］, Ｅ［ｘｘT｜ｙ；Θ´］の更新
M-step：状態系列、ノート制御パラメータ、表現制御パラメータ、上記（３４）〜（３８）式の更新
収束判定：事後確率Ｐ（Θ｜ｙ）、もしくはパラメータ値が収束していなければ、Ｅ−ｓｔｅｐへ
【００７１】
＜システム構成＞
次に、図３を参照して、本発明の実施形態による音響信号分析装置の構成を説明する。本発明の実施の形態に係る音響信号分析装置１０は、歌声を示す音響信号の基本周波数軌跡を、ノート指令入力信号と線形２次系システムによって導出されるノート制御フィルタ特性を表すインパルス応答信号との畳みこみによって得られるノート成分と、表現指令入力信号と臨界制動２次系システムによって導出される表現制御フィルタ特性を表すインパルス応答信号との畳みこみによって得られる表現成分と、微細変動成分との和で表したときの、ノート指令入力信号を表すモデルを構成するためのノート指令入力信号パラメータと、ノート制御フィルタ特性を表すモデルを構成するためのノート制御パラメータと、表現指令入力信号を表すモデルを構成するための表現指令入力信号パラメータと、表現制御フィルタ特性を表すモデルを構成するための表現制御パラメータと、前記微細変動成分を表すモデルを構成するための微細変動成分パラメータとを推定する。ここで、ノート指令入力信号は階段状の信号であり、各時刻の状態を示す状態系列パラメータと、各音符の音高からの推移を示す音高シフトパラメータとからなるノート指令入力信号パラメータを用いてモデル化され、表現指令入力信号は矩形状のパルス信号であり、各時刻の状態を示す状態系列パラメータと、各音符の各状態における表現指令の大きさを表す表現振幅パラメータとからなる表現指令入力信号パラメータを用いて、隠れマルコフモデルによって確率モデル化され、微細変動成分は、微細変動成分パラメータからなるガウス性白色雑音として確率モデル化される。図３に示すように、音響信号分析装置１０は、入力部１２と、後述するパラメータ推定処理ルーチンを実行するコンピュータ１４と、コンピュータ１４の処理により算出されたモデルパラメータを出力する出力部１６と、を備えている。
【００７２】
入力部１２により、歌声音響信号及び歌唱曲の楽譜がコンピュータ１４に入力される。
【００７３】
コンピュータ１４は、音響信号分析装置１０全体の制御を司るＣＰＵ、後述するパラメータ推定処理ルーチンのプログラム等を記憶した記憶媒体としてのＲＯＭ、ワークエリアとしてデータを一時格納するＲＡＭ、及びこれらを接続するバスを含んで構成されている。このような構成の場合には、各構成要素の機能を実現するためのプログラムをＲＯＭやＨＤＤ等の記憶媒体に記憶しておき、これをＣＰＵが実行することによって、各機能が実現されるようにする。
【００７４】
このコンピュータ１４をハードウエアとソフトウエアとに基づいて定まる機能実現手段毎に分割した機能ブロックで説明すると、図３に示すように、基本周波数抽出部１８と、パラメータ初期値データベース２０と、初期状態系列推定部２２と、信号分離部２４と、モデルパラメータ更新部２６と、パラメータ収束判定部３６と、を含んだ構成で表すことができる。
【００７５】
基本周波数抽出部１８は、入力部１２により入力される歌声音響信号から基本周波数軌跡を推定し出力する。この処理は、周知技術により実現でき、例えば、文献：A de Cheveign´e and H. Kawahara,“ YIN, a fundamental frequency estimator for speech and music,” Journal of the AcousticalSociety of America, vol.111, no.4, pp. 1917-1930, 2002 で提案される基本周波数推定法ＹＩＮを利用する。この手法は自己相関関数を使用して基本周波数を推定するが、倍ピッチエラーや半ピッチエラー、その他ノイズによる推定エラーを低減するために差分関数、正規化、放物線補間などの処理を後処理として導入した手法である。ピッチの高い音楽や歌声の基本周波数推定に効果的な手法であることが従来研究によって明らかにされている。本実施形態では、ＹＩＮを利用して、歌声音響信号から５ｍｓごとに基本周波数を推定し、基本周波数軌跡を出力する。
【００７６】
パラメータ初期値データベース２０は、入力部１２により入力された歌唱旋律の各音符を示す音高を示す音高目標値時系列、及びノート指令入力信号を表すモデルを構成するためのノート指令入力信号パラメータと、ノート制御フィルタ特性を表すモデルを構成するためのノート制御パラメータと、表現指令入力信号を表すモデルを構成するための表現指令入力信号パラメータと、表現制御フィルタ特性を表すモデルを構成するための表現制御パラメータと、微細変動成分を表すモデルを構成するための微細変動成分パラメータとからなるモデルパラメータの各々の初期値を保有するデータベースである。また、パラメータ初期値データベース２０は、その他の各種パラメータが記憶されている。まず、入力部１２により入力される歌唱曲の楽譜に応じて、そのメロディの音符数Ｉと各音符の音高{Ａn(i)}i=1Iがパラメータ初期値データベース２０に蓄積される。本実施形態の表現指令入力信号を生成するための状態数ＪはＪ＝５と固定し、パラメータ初期値データベース２０に記憶される。状態遷移確率は、
【００７７】
【数３４】
2014026055.tif 000036
に固定されパラメータ初期値データベース２０に蓄積される。事前分布におけるパラメータはμφ＝６,σφ2＝０．１,μψ＝０．６, σψ2＝０．０２,μξ＝３, σξ2＝０．１,σd2＝２５００,σB2＝１００,
【００７８】
【数３５】
2014026055.tif 000037
（１≦ｉ≦Ｉ）に固定され、パラメータ初期値データベース２０に蓄積される。音高シフトパラメータ{ｄi}i=1 I、表現指令の大きさパラメータ{Ｂe(i,j)}i=1,j=1I,Jの初期値はすべて０としてパラメータ初期値データベース２０に蓄積される。分散パラメータσn2,σe2の初期値は１００としてパラメータ初期値データベース２０に蓄積される。これらのパラメータの初期値は実験的に決定した値である。ノート制御パラメータ{φ(i),ψ(i)}i=1 Iの初期値は、上記の非特許文献３で得られた知見に基づいて、φ＝ ６．５ｒａｄ／ｍｓ，ψ＝ ０．５としてパラメータ初期値データベース２０に蓄積される。表現制御パラメータξの初期値も上記の非特許文献３で得られた知見に基づいて、ξ＝ ０．３としてパラメータ初期値データベース２０に蓄積される。
【００７９】
初期状態系列推定部２２は、基本周波数軌跡と、パラメータ初期値データベースから得られるモデルパラメータの初期値及び各音符の音高を入力として、図２のＨＭＭのビタビアライメントを行うことによって初期状態系列｛ｓk｝k=1Kを推定し、基本周波数軌跡と初期状態系列、モデルパラメータ初期値を出力する。具体的には、ｏ［ｋ］に、基本周波数軌跡のｋ番目の要素を代入して、上記（２７)式〜（２９)式に従って、初期状態系列｛ｓk｝k=1Kを推定する。
【００８０】
信号分離部２４では、ＥＭアルゴリズムに基づいて導出されるノート成分、表現成分、及び微細変動成分から構成される完全データの期待値と、完全データの自己相関とを用いて、基本周波数軌跡をノート成分と表現成分と微細変動成分に分離する。ここでは、基本周波数軌跡と現在のモデルパラメータΘ´＝｛θu,｛φ(i),ψ(i)｝i=1I,ξ,σf2｝と各音符の音高とを入力として、上記（２２）式、（２３）式に従って完全データの条件付き期待値Ｅ［ｘ｜ｙ；Θ´］，完全データの自己相関の条件付き期待値Ｅ［ｘｘT｜ｙ；Θ´］を計算し、上記（２４）式に基づいて、￣ｘn,￣ｘe, ￣ｘf，Ｒn, Ｒe, Ｒfを出力する。
【００８１】
モデルパラメータ更新部２６は、ノート制御パラメータ更新部２８と、表現制御パラメータ更新部３０と、指令入力信号パラメータ更新部３２と、微細変動成分パラメータ更新部３４とから構成されており、信号分離部２４で出力された完全データの条件付き期待値及び完全データの自己相関の条件付き期待値、並びにモデルパラメータの事前確率を用いて表される、上記（２１）式のＱ関数を最大化するように、それぞれのモデルパラメータを更新する。
【００８２】
ノート制御パラメータ更新部２８は、線形２次系システムの固有周波数に反比例するパラメータと、線形２次系システムの減衰率パラメータとからなるノート制御パラメータ｛φ(i),ψ(i)｝i=1Iの値を更新する。基本周波数軌跡、各音符の音高、信号分離部２４で出力された￣ｘn, Ｒn、そしてノート制御パラメータ以外のモデルパラメータを入力として、上記（３０）式を最大化するφ(i),ψ(i)をニュートン・ラフソン法を利用して数値的に導出し出力する。これをi = 1,...,Iに関して繰り返す。具体的には、上記（３０）式を、ノート制御パラメータの、固有周波数に反比例するパラメータと減衰率パラメータとでそれぞれ微分して得られた連立方程式を解くことにより、ノート制御パラメータの値を更新する。
【００８３】
表現制御パラメータ更新部３０は、臨界制動２次系システムの固有周波数に反比例する表現制御パラメータξの値を更新する。基本周波数軌跡、各音符の音高、信号分離部２４で出力された￣ｘe,Ｒe、そして表現制御パラメータ以外のモデルパラメータを入力として、上記（３３）式を最大化するξをニュートン・ラフソン法を利用して数値的に導出し出力する。具体的には、上記（３３）式を、表現制御パラメータの固有周波数に反比例するパラメータについて微分して得られる方程式を解くことにより、表現制御パラメータの値を更新する。
【００８４】
指令入力信号パラメータ更新部３２は、基本周波数軌跡、各音符の音高、信号分離部２４で出力された￣ｘn,￣ｘe,Ｒn,Ｒe、そして指令入力信号パラメータ以外のモデルパラメータを入力として、動的計画法（上記（２７）式、（２８）式、（２９）式）によって状態系列を更新し、上記（３４）式、（３５）式、（３６）式、（３７）式によって各指令入力信号パラメータを更新し出力する。具体的には、ノート指令入力信号パラメータと表現指令入力信号パラメータをＶｉｔｅｒｂｉ学習によって更新する。
【００８５】
微細変動成分パラメータ更新部３４は、基本周波数軌跡、信号分離部２４で出力されたＲf、そして微細変動成分パラメータ以外のモデルパラメータを入力として、上記（３８）式 に基づいて、微細変動成分パラメータを更新し出力する。具体的には、上記（３８）式を、微細変動成分パラメータで微分して得られる方程式を解くことにより、微細変動成分パラメータを更新する。
【００８６】
パラメータ収束判定部３６は、信号分離部２４で計算された￣ｘn,￣ｘe,Ｒn, Ｒe, Ｒfおよび、モデルパラメータ更新部２６によってそれぞれ更新されたモデルパラメータΘ＝｛θu,｛φ(i),ψ(i)｝i=1I,ξ,σf2｝を利用して、上記（１７）式、（１８）式の事後確率Ｐ（Θ｜ｙ）∝Ｐ（ｙ｜Θ）Ｐ（Θ）の値を計算する。更新前のモデルパラメータを用いて計算した事後確率の値と更新後のモデルパラメータを用いて計算した事後確率の値との誤差が、所定の閾値以下であれば、収束したと判定する。本実施形態ではこの誤差を実験的にε＝１．０×１０-5とする  。
【００８７】
出力部１６は、パラメータ収束判定部３６において、収束したと判定された場合にはモデルパラメータΘをすべて出力する。
【００８８】
＜音響信号分析装置の作用＞
次に、本実施の形態に係る音響信号分析装置１０の作用について説明する。まず、歌唱曲の楽譜から得られる各音符の音高及び音符数が、入力部１２により音響信号分析装置１０に入力されると、パラメータ初期値データベース２０に格納される。そして、分析対象の音響信号として、観測された当該歌唱者の歌声音響信号の時系列データが、音響信号分析装置１０に入力されると、図４に示すパラメータ推定処理ルーチンが実行される。
【００８９】
ステップＳ１０２で、基本周波数抽出部１８によって、基本周波数推定法ＹＩＮを利用して、入力された歌声音響信号から５ｍｓごとに基本周波数を推定し、基本周波数軌跡を推定する。
【００９０】
ステップＳ１０４で、上記ステップＳ１０２で推定された基本周波数軌跡と、パラメータ初期値データベース２０から得られるモデルパラメータの初期値及び各音符の音高を入力として、上記（２７）式〜（２９）式に従って、初期状態系列｛ｓk｝k=1Kを推定し、基本周波数軌跡と初期状態系列、モデルパラメータ初期値を出力する。
【００９１】
ステップＳ１０６で、上記ステップＳ１０２で推定された基本周波数軌跡と、パラメータ初期値データベース２０から得られるモデルパラメータの初期値又は後述するステップＳ１０８〜１１４で前回推定されたモデルパラメータの値と、各音符の音高とを入力として、上記（２２）式、（２４）式に従って、ノート成分、表現成分、及び微細変動成分から構成される完全データの条件付き期待値と、完全データの自己相関の条件付き期待値とを算出して、基本周波数軌跡をノート成分と表現成分と微細変動成分に分離する。
【００９２】
ステップＳ１０８で、上記ステップＳ１０６の算出結果と、パラメータ初期値データベース２０から得られるモデルパラメータの初期値又はステップＳ１０８〜１１４で前回推定されたモデルパラメータの値と、各音符の音高とを入力として、上記（３０）式をノート制御パラメータφ、ψでそれぞれ微分して得られた連立方程式を解くことにより、ノート制御パラメータの値を更新する。
【００９３】
ステップＳ１１０で、上記ステップＳ１０６の算出結果と、パラメータ初期値データベース２０から得られるモデルパラメータの初期値又はステップＳ１０８〜１１４で前回推定されたモデルパラメータの値と、各音符の音高とを入力として、上記（３３）式を表現制御パラメータζについて微分して得られる方程式を解くことにより、表現制御パラメータの値を更新する。
【００９４】
ステップＳ１１２で、上記ステップＳ１０２で推定された基本周波数軌跡と、パラメータ初期値データベース２０から得られるモデルパラメータの初期値又はステップＳ１０８〜１１４で前回推定されたモデルパラメータの値と、各音符の音高とを入力として、上記（２７）式〜（２９）式に従って、状態系列｛ｓk｝k=1Kを推定して更新する。また、上記ステップＳ１０６の算出結果と、パラメータ初期値データベース２０から得られるモデルパラメータの初期値又はステップＳ１０８〜１１４で前回推定されたモデルパラメータの値と、各音符の音高とを入力として、上記(３４)式〜（３７）式に従って、ノート指令入力信号パラメータと表現指令入力信号パラメータを更新する。
【００９５】
ステップＳ１１４で、上記ステップＳ１０６の算出結果を入力として、上記（３８）式を微細変動成分パラメータで微分して得られる方程式を解くことにより、微細変動成分パラメータを更新する。
【００９６】
ステップＳ１１６で、モデルパラメータが収束しているか判定する。モデルパラメータが収束している場合には、ステップＳ１１８で、最終的に得られたモデルパラメータの全てを、出力部１６により出力して処理を終了し、一方、モデルパラメータが収束していない場合には、上記ステップＳ１０６へ戻る。
【００９７】
＜実験結果＞
次に、本発明の効果および作用を示すため、本発明の実施形態による音響信号分析装置を用いた実験結果を以下に説明する。ここでは、上記の実施形態を実装して、基本周波数軌跡のノート成分、表現成分、微細変動成分への分解性能を評価する。
【００９８】
評価実験では、本発明が局所最適化問題を解決できているかを人工的に作成した基本周波数軌跡を用いて確認する。基本周波数軌跡（テストデータ、図５）は、上記の非特許文献３に基づいて作成した。まず、１６種類の楽曲の楽譜から抜粋した音符列を用いて、ＩＯＩ(Interonset interval) に基づいて、１６個の階段状軌跡を作成する。ここで、サンプリング周期は５ｍｓとし、オーバーシュートを畳み込むためのパラメータはζ＝０．５７０、Ω＝０．０３６３ｒａｄ／ｍｓ、プレパレーションを畳み込むためのパラメータはζ＝０．６７５、Ω＝０．０３０８ｒａｄ／ｍｓ、ビブラートを付加するためのVibrato extent とVibrato rate はそれぞれ、８３ｃｅｎｔと６．２５Ｈｚとした。各音符の音高シフトは２５ｃｅｎｔとした。最後に、変調周波数２０Ｈｚ、振幅２０ｃｅｎｔの微細変動成分を足し合わせたものを，歌声の基本周波数軌跡とし、これをテストデータとする。
図５にパラメータの推定結果例を示す。左側がテストデータであり、右側が推定結果である。左右を比較すると正解に近いパラメータが推定されることを定性的に判断できる。それぞれの二乗平均平方根誤差（Root Mean Square Error，RMSE）を計算したところ、ノート成分はＲＭＳＥ＝７．４７ｃｅｎｔ、表現成分はＲＭＳＥ＝９．５１ｃｅｎｔ であった。対数周波数上の１００ｃｅｎｔは半音に相当し、ＲＭＳＥがその１０％以下に収まっているため、パラメータ学習アルゴリズムの動作確認および局所最適化問題を解決できていることを確認できた。
【００９９】
以上説明したように、本発明の実施の形態の音響信号分析装置１０によれば、ノート制御パラメータ、表現制御パラメータ、各指令入力信号パラメータ、微細変動成分パラメータの各パラメータが収束するまで更新を繰り返すことによって、歌声を示す音響信号の基本周波数軌跡を、ノート成分と、表現成分と、微細変動成分との和で表したときの、ノート指令入力信号パラメータと、ノート制御パラメータと、表現指令入力信号パラメータと、表現制御パラメータと、微細変動成分パラメータとを推定することができ、歌唱者の表現意図による動的変動成分（ビブラートなどの表現成分）を微細変動成分と分離して、明確に特徴づけることができる。
【０１００】
なお、上記の実施の形態では、モデルパラメータが収束したか否かを判定する方法として、更新前のモデルパラメータを用いて計算した事後確率の値と更新後のモデルパラメータを用いて計算した事後確率の値との誤差が、所定の閾値以下であれば、収束したと判定する場合を例に説明したが、これに限定されるものではない。例えば、収束したか否かを判定する方法としては、事後確率を用いる方法以外に、モデルパラメータ各々の値を更新前と更新後とで比較しても良いし、予め定めた繰り返し回数に到達したか否かで判定を行っても良い。本実施形態ではモデルパラメータ各々の値を更新前と更新後とで比較する場合、この誤差がε＝１．０×１０-5であれば良好な結果であることを実験的に確認した。また、予め定めた繰り返し回数に到達したか否かで判定をする場合、１０００回の繰り返し回数が必要であることも実験的に確認している。
【０１０１】
また、本発明は、上記実施の形態に限定されるものではなく、この発明の要旨を逸脱しない範囲内で様々な変形や応用が可能である。
【０１０２】
また、上述の音響信号分析装置１０は、内部にコンピュータシステムを有しているが、「コンピュータシステム」は、ＷＷＷシステムを利用している場合であれば、ホームページ提供環境（あるいは表示環境）も含むものとする。
【０１０３】
また、本願明細書中において、プログラムが予めインストールされている実施形態として説明したが、当該プログラムを、コンピュータ読み取り可能な記録媒体に格納して提供することも可能であるし、ネットワークを介して提供することも可能である。また、本実施の形態の音響信号分析装置１０の各部をハードウエアにより構成してもよい。また、パラメータ初期値が記憶されるデータベースとしては、ハードディスク装置やファイルサーバ等に例示される記憶手段によって実現可能であり、音響信号分析装置１０内部にデータベースを設けても良いし、外部装置に設けてもよい。
【符号の説明】
【０１０４】
１０   音響信号分析装置
１２   入力部
１４   コンピュータ
１６   出力部
１８   基本周波数抽出部
２０   パラメータ初期値データベース
２２   初期状態系列推定部
２４   信号分離部
２６   モデルパラメータ更新部
２８   ノート制御パラメータ更新部
３０   表現制御パラメータ更新部
３２   指令入力信号パラメータ更新部
３４   微細変動成分パラメータ更新部
３６   パラメータ収束判定部
【図１】
2014026055.tif 000038
【図２】
2014026055.tif 000039
【図３】
2014026055.tif 000040
【図４】
2014026055.tif 000041
【図５】
2014026055.tif 000042
【図６】
2014026055.tif 000043
【図７】
2014026055.tif 000044