(57)【要約】【課題】楽器音一音単位に信号を分解すると共に、時間軸方向の分割と、楽器音が同時出現する場合とを考慮した構造解析を実現することができるようにする。【解決手段】木構造サンプリング部３１によって、音響信号の時系列データを頂点ノードとして、時間軸方向に音響信号を分割する次元と同時に発音する楽器音に分割する次元とにおいて分割した木構造の各分割位置のシンボル及び発音時刻を、サンプリングにより決定する。パラメータサンプリング部３２によって、Ｎ個の楽器音ｎの各々に対する時間周波数成分Ｃω,t,n、音量Ｖn、発音時刻に関するパラメータτn、各シンボルkに対する音の長さＬk、シンボルｋを持つ親ノードからシンボルｉを持つ子ノード及びシンボルｊを持つ子ノードに分岐させる重みを表わすパラメータＺi,j(k)、及び各シンボルｋの重みを表わすパラメータβkの各値を、サンプリングによって決定する。【選択図】図２
2013195575.tif 000002
【特許請求の範囲】
【請求項１】
複数の楽器音が混在する音響信号の時系列データを入力として、観測時間周波数成分Ｙω,t（ωは周波数、ｔは時刻のインデックスである。）を要素にもつ二次元配列Ｙ＾を出力する時間周波数分解手段と、
Ｎ個の楽器音ｎの各々に対する時間周波数成分Ｃω,t,nを要素にもつ三次元配列Ｃ＾、前記Ｎ個の楽器音ｎの各々に対する音量Ｖn、シンボルｋが各ノードに付与され、かつ、前記音響信号の時系列データを頂点ノードとして、時間軸方向に前記音響信号を分割する次元と同時に発音する楽器音に分割する次元とにおいて、末端ノードが前記楽器音の一音に対応するように分割した木構造の各分割位置における各子ノードの発音時刻ｔChildleft、ｔChildright、前記Ｎ個の楽器音ｎの各々に対する発音時刻に関するパラメータτn、各シンボルkに対する音の長さＬk、シンボルｋを持つ親ノードからシンボルｉを持つ子ノード及びシンボルｊを持つ子ノードに分岐させる重みを表わすパラメータＺi,j(k)、及び各シンボルｋの重みを表わすパラメータβkの各々の初期値を設定するパラメータ初期値設定手段と、
前記子ノードの発音時刻ｔChildleft、ｔChildright、前記音の長さＬk、及び前記パラメータＺi,j(k)に基づいて求められる、前記木構造の各分割位置ｍにおける各子ノードのシンボルChildleft、Childright、発音時刻ｔChildleft、ｔChildrightの確率分布に従って、前記木構造の各分割位置ｓにおける各子ノードのシンボルChildleft、Childright、発音時刻ｔChildleft、ｔChildrightをサンプリングにより決定することにより前記木構造を生成する木構造サンプリング手段と、
前記木構造サンプリング手段によって生成された前記木構造に基づいて前記パラメータＺi,j(k)の各値をサンプリングすると共に、前記三次元配列Ｃ＾、前記音量Ｖn、前記発音時刻に関するパラメータτn、前記音の長さＬk、及び前記パラメータβkの各々の値をサンプリングして、前記時間周波数分解手段によって出力された前記二次元配列Ｙ＾、前記木構造、前記三次元配列Ｃ＾、前記音量Ｖn、前記発音時刻に関するパラメータτn、前記音の長さＬk、前記パラメータＺi,j(k)、及び前記パラメータβkに基づいて求められる対数事後分布に従って、前記パラメータＺi,j(k)、前記三次元配列Ｃ＾、前記音量Ｖn、前記発音時刻に関するパラメータτn、前記音の長さＬk、及び前記パラメータβkの各々の値を更新するパラメータサンプリング手段と、
予め定められた終了条件を満たすまで、前記木構造サンプリング手段による生成及び前記パラメータサンプリング手段による更新を反復して行い、反復する毎に前記パラメータサンプリング手段によって得られた前記三次元配列Ｃ＾の前記Ｎ個の楽器音ｎの各々に対する時間周波数成分Ｃω,t,nを出力する終了条件判定手段と、
を含む音響信号分析装置。
【請求項２】
前記パラメータサンプリング手段は、さらに前記対数事後分布に基づいて前記更新された前前記パラメータＺi,j(k)、記三次元配列Ｃ＾、前記音量Ｖn、前記発音時刻に関するパラメータτn、前記音の長さＬk、及び前記パラメータβkに対する確率を出力し、
前記終了条件判定手段によって出力された、反復する毎に前記パラメータサンプリング手段によって得られた前記三次元配列Ｃ＾のうち、最も高い前記確率に対応する前記三次元配列Ｃ＾を選択する出力結果選択手段を更に含む請求項１記載の音響信号分析装置。
【請求項３】
時間周波数分析手段、パラメータ初期値設定手段、木構造サンプリング手段、パラメータサンプリング手段、及び終了条件判定手段を含む音響信号分析装置における音響信号分析方法であって、
前記時間周波数分解手段によって、複数の楽器音が混在する音響信号の時系列データを入力として、観測時間周波数成分Ｙω,t（ωは周波数、ｔは時刻のインデックスである。）を要素にもつ二次元配列Ｙ＾を出力し、
前記パラメータ初期値設定手段によって、Ｎ個の楽器音ｎの各々に対する時間周波数成分Ｃω,t,nを要素にもつ三次元配列Ｃ＾、前記Ｎ個の楽器音ｎの各々に対する音量Ｖn、シンボルｋが各ノードに付与され、かつ、前記音響信号の時系列データを頂点ノードとして、時間軸方向に前記音響信号を分割する次元と同時に発音する楽器音に分割する次元とにおいて、末端ノードが前記楽器音の一音に対応するように分割した木構造の各分割位置における各子ノードの発音時刻ｔChildleft、ｔChildright、前記Ｎ個の楽器音ｎの各々に対する発音時刻に関するパラメータτn、各シンボルkに対する音の長さＬk、シンボルｋを持つ親ノードからシンボルｉを持つ子ノード及びシンボルｊを持つ子ノードに分岐させる重みを表わすパラメータＺi,j(k)、及び各シンボルｋの重みを表わすパラメータβkの各々の初期値を設定し、
前記木構造サンプリング手段によって、前記子ノードの発音時刻ｔChild、前記音の長さＬk、及び前記パラメータＺi,j(k)に基づいて求められる、前記木構造の各分割位置ｍにおける各子ノードのシンボルChildleft、Childright、発音時刻ｔChildleft、ｔChildrightの確率分布に従って、前記木構造の各分割位置ｓにおける各子ノードのシンボルChildleft、Childright、発音時刻ｔChildleft、ｔChildrightをサンプリングにより決定することにより前記木構造を生成し、
前記パラメータサンプリング手段によって、前記木構造サンプリング手段によって生成された前記木構造に基づいて前記パラメータＺi,j(k)の各値をサンプリングすると共に、前記三次元配列Ｃ＾、前記音量Ｖn、前記発音時刻に関するパラメータτn、前記音の長さＬk、及び前記パラメータβkの各々の値をサンプリングして、前記時間周波数分解手段によって出力された前記二次元配列Ｙ＾、前記木構造、前記三次元配列Ｃ＾、前記音量Ｖn、前記発音時刻に関するパラメータτn、前記音の長さＬk、前記パラメータＺi,j(k)、及び前記パラメータβkに基づいて求められる対数事後分布に従って、前記パラメータＺi,j(k)、前記三次元配列Ｃ＾、前記音量Ｖn、前記発音時刻に関するパラメータτn、前記音の長さＬk、及び前記パラメータβkの各々の値を更新し、
前記終了条件判定手段によって、予め定められた終了条件を満たすまで、前記木構造サンプリング手段による生成及び前記パラメータサンプリング手段による更新を反復して行い、反復する毎に前記パラメータサンプリング手段によって得られた前記三次元配列Ｃ＾の前記Ｎ個の楽器音ｎの各々に対する時間周波数成分Ｃω,t,nを出力する
音響信号分析方法。
【請求項４】
コンピュータを、請求項１又は２記載の音響信号分析装置の各手段として機能させるためのプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、音響信号分析装置、方法、及びプログラムに係り、特に、音響信号の時系列から、信号パラメータを分析する音響信号分析装置、方法、及びプログラムに関する。
【背景技術】
【０００２】
従来、楽器音の重畳した信号から各楽器音ごとに信号を分解する手法に非負値行列因子分解(Non-negative matrix factorization, NMF)による方法が知られている（例えば、非特許文献１）。NMFは、信号の時間周波数表現であるスペクトログラムY＾=(Yω,t)Ω×T∈R≧0（ただし、ω=1,...,Ωは周波数のインデックス、t=1,...,Tは時間のインデックスを表す）の近似モデルとしてN個のスペクトルH＾=(Hω,n)Ω×Nと、それらのゲインU＾=(Un,t)N×Tによって、
【０００３】
【数１】
2013195575.tif 000003
【０００４】
となるH＾、U＾を求める問題に帰着される。観測スペクトログラムと近似モデルの間の誤差を測る尺度はさまざまなものが用いられているが、ここではその一例として、一般化カルバックライブラーダイバージェンスを用いる。NMFはコンポーネントと呼ばれる隠れた単一楽器音スペクトログラムCn＾を用いて
【０００５】
【数２】
2013195575.tif 000004
【０００６】
と表せることが知られている。NMFを用いる上でNの値の設定は大きな問題であり、事前に指定することなく観測データから自動的に推定される枠組みが求められてきた。そこでNMFは下記のようにGamma processを事前分布にしたBayesモデルとして取り扱われることがある。
【０００７】
【数３】
2013195575.tif 000005
【０００８】
ただし、η、λは指定するパラメータである。これによって、Nは十分に大きな値に設定しておけばNの設定が分解の能力に大きな影響を与えないことが知られている。Nを正の無限大に近づけるに従ってモデルの近似誤差が小さくなるので、計算機環境、計算時間の許す限り大きな値を用いるのがよい。なお、記号に付された「＾」は、当該記号が行列または多次元配列またはベクトルであることを表わしている。
【０００９】
また、従来、系列データの解析に確率文脈自由文法(Probabilistic context-free grammar, PCFG)がよく用いられている。PCFGは文脈自由文法（シンボルと分岐規則A→B C or A→w,ただしA,B,Cは非終端記号でwは終端記号）と分岐規則の確率分布によって定義され、木構造への確率モデルと見なすことが出来る。従来からシンボルの数や分岐規則の数に制約を置かず観測データからそれらも推定する枠組みとしてinfinite PCFGが提案されており良く用いられている（例えば、非特許文献２）。infinite PCFGは階層ディリクレプロセス(Hierarchical Dirichlet process,HDP)を用いて次のような手順で構成される。
【００１０】
まず、シンボルの領域に対して次のようなディリクレプロセス(Dirichlet process, DP)から生成された確率測度を考える。
【００１１】
【数４】
2013195575.tif 000006
【００１２】
ただし、Iはシンボル領域への基底測度でγはDPの集中度パラメータであり、A〜BはAが確率分布Bから生成されたことを表す。基底測度は直感的にはシンボルの出現しやすさを大まかに表す確率分布であると考えればよい。集中度γは、DPから生成された確率測度と基底測度をどの程度似せるかを表すパラメータであると見なすことが出来る。ディリクレ過程から生成された確率測度G0は、基底測度に似た傾向の離散分布をシンボル領域に構成する。G0は可算無限個のシンボルの候補（アトム）φk(k=1,2,...)とそれらがどの程度出現しやすいかを表す重みβ（合計して1になる無限次元のベクトル）であると見なすことが出来る。k番目のインデックスを持つシンボルからの分岐規則に関する確率分布はββTを基底測度とするG'kから生成されたとみなすことが出来る。
【００１３】
【数５】
2013195575.tif 000007
【００１４】
ただし、G'0=Σi,jβiβjδ(φi,φj)．である。φiはシンボル領域のアトムを表している。各G'kは直感的にはk番目のシンボルから分岐によって生成されるシンボルの組に対する無限次元多項分布と考えることが出来る。全てのG'k(k=1,2,...)はアトムを共有するため、HDPによる事前分布を設定したことによって可算無限個のシンボルの中で分岐していく木構造への確率分布を構成することが出来る。
【先行技術文献】
【非特許文献】
【００１５】
【非特許文献１】M. Hoffman, D. Blei and P. Cook, \Bayesian nonparametric matrix factorization for recorded music," ICML, 2010.
【非特許文献２】P. Liang, S. Petrov, M. I. Jordan and D. Klein, \The in_nite PCFG using hierarchical Dirichlet processes," EMNLP, 2007.
【発明の概要】
【発明が解決しようとする課題】
【００１６】
従来の信号分解では各楽器音ごとに信号を分解する際に、その音符情報を得ることは出来なかった。例えば、ピアノのCの音3回、ヴァイオリンのDの音2回が重畳した信号を考える。この時、NMFはこの信号をピアノのCとヴァイオリンのDとの２つの信号及び音量に分解する働きをするが、それらがピアノのＣの音３回、ヴァイオリンのＤの音２回演奏されていることを特定することは出来ていなかった。しかし各楽器音の出現位置を検出するためには、各楽器音がいつ何回演奏されたのかを推定することが不可欠になる。
【００１７】
また、従来の構造解析は一次元の階層構造しか扱うことが出来なかった。複数の楽器音が重畳する信号においては、音楽の和音のように同時に出現する楽器音が存在する。
【００１８】
本発明は、上記の事情を考慮してなされたもので、楽器音一音単位に信号を分解すると共に、時間軸方向の分割と、楽器音が同時出現する場合とを考慮した構造解析を実現することができる音響信号分析装置、方法、及びプログラムを提供することを目的とする。
【課題を解決するための手段】
【００１９】
上記の目的を達成するために本発明に係る音響信号分析装置は、複数の楽器音が混在する音響信号の時系列データを入力として、観測時間周波数成分Ｙω,t（ωは周波数、ｔは時刻のインデックスである。）を要素にもつ二次元配列Ｙ＾を出力する時間周波数分解手段と、Ｎ個の楽器音ｎの各々に対する時間周波数成分Ｃω,t,nを要素にもつ三次元配列Ｃ＾、前記Ｎ個の楽器音ｎの各々に対する音量Ｖn、シンボルｋが各ノードに付与され、かつ、前記音響信号の時系列データを頂点ノードとして、時間軸方向に前記音響信号を分割する次元と同時に発音する楽器音に分割する次元とにおいて、末端ノードが前記楽器音の一音に対応するように分割した木構造の各分割位置における各子ノードの発音時刻ｔChildleft、ｔChildright、前記Ｎ個の楽器音ｎの各々に対する発音時刻に関するパラメータτn、各シンボルkに対する音の長さＬk、シンボルｋを持つ親ノードからシンボルｉを持つ子ノード及びシンボルｊを持つ子ノードに分岐させる重みを表わすパラメータＺi,j(k)、及び各シンボルｋの重みを表わすパラメータβkの各々の初期値を設定するパラメータ初期値設定手段と、前記子ノードの発音時刻ｔChild、前記音の長さＬk、及び前記パラメータＺi,j(k)に基づいて求められる、前記木構造の各分割位置ｍにおける各子ノードのシンボルChildleft、Childright、発音時刻ｔChildleft、ｔChildrightの確率分布に従って、前記木構造の各分割位置ｓにおける各子ノードのシンボルChildleft、Childright、発音時刻ｔChildleft、ｔChildrightをサンプリングにより決定することにより前記木構造を生成する木構造サンプリング手段と、前記木構造サンプリング手段によって生成された前記木構造に基づいて前記パラメータＺi,j(k)の各値をサンプリングすると共に、前記三次元配列Ｃ＾、前記音量Ｖn、前記発音時刻に関するパラメータτn、前記音の長さＬk、及び前記パラメータβkの各々の値をサンプリングして、前記時間周波数分解手段によって出力された前記二次元配列Ｙ＾、前記木構造、前記三次元配列Ｃ＾、前記音量Ｖn、前記発音時刻に関するパラメータτn、前記音の長さＬk、前記パラメータＺi,j(k)、及び前記パラメータβkに基づいて求められる対数事後分布に従って、前記パラメータＺi,j(k)、前記三次元配列Ｃ＾、前記音量Ｖn、前記発音時刻に関するパラメータτn、前記音の長さＬk、及び前記パラメータβkの各々の値を更新するパラメータサンプリング手段と、予め定められた終了条件を満たすまで、前記木構造サンプリング手段による生成及び前記パラメータサンプリング手段による更新を反復して行い、反復する毎に前記パラメータサンプリング手段によって得られた前記三次元配列Ｃ＾の前記Ｎ個の楽器音ｎの各々に対する時間周波数成分Ｃω,t,nを出力する終了条件判定手段と、を含んで構成されている。
【００２０】
本発明に係る音響信号分析方法は、時間周波数分析手段、パラメータ初期値設定手段、木構造サンプリング手段、パラメータサンプリング手段、及び終了条件判定手段を含む音響信号分析装置における音響信号分析方法であって、前記時間周波数分解手段によって、複数の楽器音が混在する音響信号の時系列データを入力として、観測時間周波数成分Ｙω,t（ωは周波数、ｔは時刻のインデックスである。）を要素にもつ二次元配列Ｙ＾を出力し、前記パラメータ初期値設定手段によって、Ｎ個の楽器音ｎの各々に対する時間周波数成分Ｃω,t,nを要素にもつ三次元配列Ｃ＾、前記Ｎ個の楽器音ｎの各々に対する音量Ｖn、シンボルｋが各ノードに付与され、かつ、前記音響信号の時系列データを頂点ノードとして、時間軸方向に前記音響信号を分割する次元と同時に発音する楽器音に分割する次元とにおいて、末端ノードが前記楽器音の一音に対応するように分割した木構造の各分割位置における各子ノードの発音時刻ｔChildleft、ｔChildright、前記Ｎ個の楽器音ｎの各々に対する発音時刻に関するパラメータτn、各シンボルkに対する音の長さＬk、シンボルｋを持つ親ノードからシンボルｉを持つ子ノード及びシンボルｊを持つ子ノードに分岐させる重みを表わすパラメータＺi,j(k)、及び各シンボルｋの重みを表わすパラメータβkの各々の初期値を設定し、前記木構造サンプリング手段によって、前記子ノードの発音時刻ｔChild、前記音の長さＬk、及び前記パラメータＺi,j(k)に基づいて求められる、前記木構造の各分割位置ｍにおける各子ノードのシンボルChildleft、Childright、発音時刻ｔChildleft、ｔChildrightの確率分布に従って、前記木構造の各分割位置ｓにおける各子ノードのシンボルChildleft、Childright、発音時刻ｔChildleft、ｔChildrightをサンプリングにより決定することにより前記木構造を生成し、前記パラメータサンプリング手段によって、前記木構造サンプリング手段によって生成された前記木構造に基づいて前記パラメータＺi,j(k)の各値をサンプリングすると共に、前記三次元配列Ｃ＾、前記音量Ｖn、前記発音時刻に関するパラメータτn、前記音の長さＬk、及び前記パラメータβkの各々の値をサンプリングして、前記時間周波数分解手段によって出力された前記二次元配列Ｙ＾、前記木構造、前記三次元配列Ｃ＾、前記音量Ｖn、前記発音時刻に関するパラメータτn、前記音の長さＬk、前記パラメータＺi,j(k)、及び前記パラメータβkに基づいて求められる対数事後分布に従って、前記パラメータＺi,j(k)、前記三次元配列Ｃ＾、前記音量Ｖn、前記発音時刻に関するパラメータτn、前記音の長さＬk、及び前記パラメータβkの各々の値を更新し、前記終了条件判定手段によって、予め定められた終了条件を満たすまで、前記木構造サンプリング手段による生成及び前記パラメータサンプリング手段による更新を反復して行い、反復する毎に前記パラメータサンプリング手段によって得られた前記三次元配列Ｃ＾の前記Ｎ個の楽器音ｎの各々に対する時間周波数成分Ｃω,t,nを出力する。
【００２１】
本発明に係るプログラムは、上記の音響信号分析装置の各手段としてコンピュータを機能させるためのプログラムである。
【発明の効果】
【００２２】
以上説明したように、本発明の音響信号分析装置、方法、及びプログラムによれば、音響信号の時系列データを頂点ノードとして、時間軸方向に前記音響信号を分割する次元と同時に発音する楽器音に分割する次元とにおいて、末端ノードが楽器音の一音に対応するように分割した木構造のサンプリングと、Ｎ個の楽器音ｎの各々に対する時間周波数成分Ｃω,t,n、音量Ｖn、発音時刻に関するパラメータτn、各シンボルkに対する音の長さＬk、シンボルｋを持つ親ノードからシンボルｉを持つ子ノード及びシンボルｊを持つ子ノードに分岐させる重みを表わすパラメータＺi,j(k)、及び各シンボルｋの重みを表わすパラメータβkのサンプリングとを反復して行うことにより、楽器音一音単位に信号を分解すると共に、時間軸方向の分割と、楽器音が同時出現する場合とを考慮した構造解析を実現することができる、という効果が得られる。
【図面の簡単な説明】
【００２３】
【図１】二次元infinite PCFGを説明するための図である。
【図２】本発明の第１の実施の形態に係る音響信号分析装置の構成を示す概略図である。
【図３】本発明の第１の実施の形態に係る音響信号分析装置における音響信号分析処理ルーチンの内容を示すフローチャートである。
【図４】本発明の第２の実施の形態に係る音響信号分析装置の構成を示す概略図である。
【図５】（Ａ）楽器音の発音位置の正解としてのピアノロールを示す図、及び（Ｂ）入力された音響信号から推定された音高及び音量情報を示す図である。
【図６】（Ａ）入力された音響信号からサンプリングされた木構造を示す図、及び（Ｂ）入力された音響信号からサンプリングされた木構造を示す図である。
【発明を実施するための形態】
【００２４】
以下、図面を参照して本発明の実施の形態を詳細に説明する。本発明で提案する手法では、楽器音一音単位からのスペクトログラムの生成モデルを構成し、それらの出現位置に対して、二次元に拡張したinfinite PCFGの事前分布を設定することによって、楽器音一音単位の信号の分解と、それらの出現位置に関する構造解析を実現する。
【００２５】
＜スペクトログラムの生成モデル＞
まず、本発明で提案する楽器音一音単位からのスペクトログラムの生成モデルを説明する。各楽器音一音は、定常なスペクトル、音量、音量変化パターン、及び発音時刻の4つの要素から構成されていると考える。
【００２６】
このときn番目の楽器音はdn番目のスペクトルhdnを一様分布から選び、同様にen番目の音量変化パターンOenを一様分布から選んだとみなす。音量をVn、発音時刻に関するパラメータをτnとすると、n番目の楽器音を表すスペクトログラムは
【００２７】
【数６】
2013195575.tif 000008
のように表せる。Poisson(a)はパラメータaのポアソン分布を表している。音量Vn、発音時刻に関するパラメータτnは、楽器音一音ずつ与えられる。
【００２８】
次に各楽器音の出現位置の構造を司る二次元infinite PCFGについて、図１を用いて説明する。これは各分岐の際に、どの次元に関して分岐を行うのかを毎回選択することによって実現することができ、音符情報に相当するシンボルｋが各ノードに付与され、かつ、音響信号の時系列データを頂点ノードとして、時間軸方向に分割する次元と同時に発音する楽器音に分割する次元とにおいて、末端ノードが楽器音一音に対応するように分割した木構造を生成する。まずは時間軸方向に分割する次元での分岐規則に関する確率モデルについて説明する。従来のinfinite PCFGと同様、シンボル領域のアトムを構成する際にDPを用いるが、このときに時間軸上での分割を扱うために、各アトムに長さ（実数）の概念を付加したい。これは
【００２９】
【数７】
2013195575.tif 000009
のようにシンボル領域への基底測度Iと実数領域への基底測度Fの直積を基底測度とするDPによって実現される。i番目のシンボルには隠れた実数の音の長さLi(i=1,2,...)が対応しているとし、従来のinfinite PCFG同様
【００３０】
【数８】
2013195575.tif 000010
としてG'k(k=1,2,...)を構成した後、親ノードのシンボルの長さが子ノードのシンボルの長さの合計に近いほど重みが高くなる係数(w i,j (k))K×KKを用いて
【００３１】
【数９】
2013195575.tif 000011
とすることで各Gkは時間軸を分割していく規則ほど使われやすくなるような分岐規則の確率モデルを構成することが出来る。ただし、∝は比例を表す。
【００３２】
次に、同時に出現する楽器音に分割する次元の分岐規則に対する確率モデルを説明する。どちらの次元への分岐規則を採用するかはm番目の分岐についてbm(bm〜Bernoulli(aB)、ただしaBは実数のパラメータ)を用い、bm=1のときはm番目の分岐は、時間軸方向に分割する次元でGkにより行われ、bm=0のときはm番目の分岐は、同時に出現する楽器音に分割する次元で、親ノードのシンボルと同じシンボルを持つ２つの子ノードを同じ時間軸上の位置に確率1で配置する。
【００３３】
このようなモデルの構成法は一通りではないが、例えば
【００３４】
【数１０】
2013195575.tif 000012
のように構成することが出来る。ただし、δは、デルタ関数を表わす。Ｚi,j(k)は、シンボルｋを持つ親ノードから各々シンボルｉ，ｊを持つ子ノードに分岐する重みを表わすパラメータである。Kはシンボル数に関して十分に大きな打ち切りの数(Kが正の無限大に近づくに従ってモデルの近似誤差が小さくなるため計算機環境の許す限り大きな値を設定するのが良い)、ChildleftとChildleftはそれぞれ子ノードのシンボルのインデックス、tleftとtrightは子ノードの発音時刻、tparentは親ノードの発音時刻である。Dirichlet、Gamma、Normalはそれぞれディリクレ分布、ガンマ分布(一つ目のパラメータが形状母数、二つ目のパラメータが逆尺度母数)、正規分布（一つ目のパラメータが平均、二つ目のパラメータが分散）を表している。確率変数の推定に関しては、tn(n=1，2,...)を連続値として考える場合は困難であるが、量子化し離散値にするとマルコフ連鎖モンテカルロ法によって行うことが出来る。具体的な手順の一例を以下で説明する。信号を短時間フーリエ変換によって変換した振幅スペクトログラムY＾を入力として、正の実数α、γ、ρと自然数K、Nを指定し、以下の（１３）に示すような対数事後分布（ただし定数項は除いた）に従う木構造の各分割位置ｍの子ノードのシンボルのインデックスの集合Childleft=(Childleft(m))M、Childright=(Childright(m))M、及び木構造の各分割位置ｍにおける子ノードの発音時刻の集合t＾={(tChildleft(m))M,(tChildright(m))M}と、Ｎ個の楽器音ｎの各々における時間周波数成分Ｃω,t,nを要素にもつ三次元配列C＾=(Cω,t,n)ω×T×Nと、Ｎ個の楽器音ｎの各々における音量V＾=(Vn)Nと、Ｎ個の楽器音ｎの各々における発音時刻に関するパラメータτ＾=(τn)Nと、各シンボルｋにおける音の長さL＾=(Lk)Kと、 重みを表わすパラメータZ＾=(Z(k)i,j)KK,K, 各シンボルｋの重みを表わすパラメータβ＾=(βk)Kのサンプルを、マルコフ連鎖モンテカルロ法によって求める。
【００３５】
【数１１】
2013195575.tif 000013
【００３６】
ただし、Parent(m),Childleft(m),Childright(m)はそれぞれ木構造のm番目の分割位置における親ノードのシンボルインデックス、左側の子ノードのシンボルインデックス、右側の子ノードのシンボルインデックスを表し、I[A=B]はA=Bのときに1となりそれ以外で0となる関数を表している。
【００３７】
＜システム構成＞
次に、音響信号の信号パラメータを分析して出力する音響信号分析装置に、本発明を適用した場合を例にして、本発明の実施の形態を説明する。
【００３８】
図２に示すように、第１の実施の形態に係る音響信号分析装置は、ＣＰＵと、ＲＡＭと、後述する音響信号分析処理ルーチンを実行するためのプログラムを記憶したＲＯＭとを備えたコンピュータで構成され、機能的には次に示すように構成されている。
【００３９】
音響信号分析装置は、時間周波数解析部１と、パラメータ初期値設定部２と、信号分析部３と、記憶部４と、出力部５とを備えている。また、信号分析部３は、木構造サンプリング部３１と、パラメータサンプリング部３２と、終了判定部３３とを備えている。
【００４０】
時間周波数解析部１は、時系列信号としての観測された、複数の楽器音が混在する音響信号x={x1,x2,...xR}（ただしRは入力信号長）を入力として、時間周波数成分（観測時間周波数成分）Ｙω,t（ω＝１,・・・,Ω,ｔ＝１,・・・,Ｔは、それぞれ周波数、時刻に対応するインデックスを示す。）を各（ω，ｔ）の要素にもつ二次元配列Ｙ＾を計算し、信号分析部３に出力する。また、計算した時間周波数成分Ｙω,tを、記憶部４に記憶しておく。より詳細には、時間周波数解析部１は、時系列信号ｙを入力として、短時間フーリエ変換（Short-Time Fourier Transform；ＳＴＦＴ）を用いて時間周波数解析を行うことにより時間周波数成分Ｙω,tを計算し、時間周波数成分Ｙω，ｔを格納した行列（振幅スペクトログラム）Ｙ＾＝（Ｙω，ｔ）Ω×Ｔを出力する。例えば、サンプリングレート16000Hzの信号xを窓幅64ms、窓シフト長32ms、窓関数としてハニング窓を用いた短時間フーリエ変換をして、Ｙ＾を得る。なお、時間周波数成分Ｙω,tは、ウェーブレット変換を用いて計算してもよい。また、複数の楽器音は、複数種類の楽器の音から構成されもよく、１種類の楽器の音から構成されてもよい。
【００４１】
パラメータ初期値設定部２は、あらかじめ決めた正の実数α、γ、ρと自然数K、Nを設定すると共に、入力時間周波数表現Y＾からパラメータC＾=(Cω,t,n)Ω×T×N, V＾=(Vn)N, t＾={(tChildleft(m))M,(tChildright(m))M}, τ＾=(τn)N,L＾=(Lk)K,Z＾=(Z(k)i,j)KK,K,β＾=(βk)Kの各値について、全て非負値の乱数により初期値を定める。
【００４２】
出力部５は、信号分析部３で求めた各パラメータC＾=(Cω,t,n)Ω×T×Nを出力する。
【００４３】
次に、信号分析部３の具体処理について説明する。
【００４４】
信号分析部３は、終了判定部３３で終了したと判定されるまで、木構造サンプリング部３１、及びパラメータサンプリング部３２の処理を反復して行う。以下では、一反復の実行における各部の処理について説明をする。
【００４５】
木構造サンプリング部３１は、まず、記憶部４に記憶されている、初期値が設定された、あるいは前回の反復の実行で得られたパラメータt＾={(tChildleft(m))M,(tChildright(m))M}, L＾=(Lk)K, Z＾=(Z(k)i,j)KK,Kを用いて、以下の（１４）式に従って、四次元の配列(pn,n',k,t)N×N×K×Tの各値を再帰的にそれぞれ計算する。
【００４６】
【数１２】
2013195575.tif 000014
ただし、sは補助変数であり木の分割位置のインデックスを表している。
【００４７】
また、木構造サンプリング部３１は、木構造のｓ番目の分割位置について、以下の（１５）式に示す多項分布に基づいて求められる、kl,kr,tl,trの全ての組み合わせの各々の確率に従って、各親ノードに対する左右の子ノードのシンボル及び発音時刻をサンプリングにより決定して、木構造を生成する。
【００４８】
【数１３】
2013195575.tif 000015
【００４９】
ただし、p(k→(kl,kr),t→(tl,tr))は、シンボルがkである親ノードからkl,krのシンボルを各々持つ２つの子ノードを生成し、かつ親ノードおよび各子ノードの発音時刻がt,tl,trとなる確率を表している。sは補助変数であり木の分割位置のインデックスを表している。m番目の分割位置に対してサンプリングにより生成される子ノードのシンボルをChildleft(m)=kl,Childright(m)=krとし、それぞれの発音時刻をtChildleft(m)=tl, tChildright(m)=trとする。
【００５０】
このように、木構造サンプリング部３１は、木構造の各分割位置ｍについて、Childleft(m)、Childright(m)、tChildleft(m)、tChildright(m)をサンプリングにより更新し、子ノードのシンボルの集合Child={(Childleft(m))M,(Childright(m))M}、子ノードの発音時刻の集合t＾={(tChildleft(m))M,(tChildright(m))M}を記憶部４に記憶しておく。
【００５１】
また、パラメータサンプリング部３２は、記憶部４に記憶されている、二次元配列Ｙ＾、前回の反復の実行で得られたパラメータC＾=(Cω,t,n)Ω×T×N, V＾=(Vn)N, τ＾=(τn)N, L＾=(Lk)K, Z＾=(Z(k)i,j)KK,K, β＾=(βk)Kと、木構造サンプリング部３１により直前の実行で求められた木構造のChild＾={(Childleft(m))M,(Childright(m))M}、t＾={(tChildleft(m))M,(tChildright(m))M}とを用いて、以下の手順に従い、C＾=(Cω,t,n)Ω×T×N,V＾=(Vn)N,τ＾=(τn)N,L＾=(Lk)K,Z＾=(Z(k)i,j)KK,K,β＾=(βk)Kの各値をそれぞれ更新する。
【００５２】
まず、Z(k)i,jの各値ついては、以下の（１７）式に示すガンマ分布から乱数によりサンプリングする。
【００５３】
【数１４】
2013195575.tif 000016
【００５４】
ただし、N i,j(k)はk番目のシンボルを持つ親ノードからi,j番目のシンボルを持つ子ノードへ分岐した回数を表しており、木構造サンプリング部３１により今回の反復で求められた木構造に基づいて求められる。このサンプリングはw(k)i,jが0に近いときに計算機上での取扱いが難しい。その場合は同じ平均値を持つ正規分布からサンプリングする。
【００５５】
そして、サンプリングにより得られた値の採択可否を、min(1, exp J'/exp J)で表される確率に従って決める。ただし、J は一反復前の値を用いたときの対数事後分布に基づく確率であり、二次元配列Ｙ＾と、一反復前に更新された値Z i,j(k)と、木構造サンプリング部３１により直前に求められた木構造と、一反復前に更新された他のパラメータの値とに基づいて、上記（１３）式に従って算出する。J'は新たにサンプリングした値を用いたときの対数事後分布に基づく確率を表しており、二次元配列Ｙ＾と、サンプリングにより求められた値Z i,j(k)と、木構造サンプリング部３１により一反復前に求められた木構造と、一反復前に更新された他のパラメータの値とに基づいて、上記（１３）式に従って算出する。棄却された場合は一反復前のサンプルをそのまま用いる。
【００５６】
他のパラメータ（Cω,t,n, Vn, τn, Lk, β）の各々に関しては、一反復前の値を平均値とする正規分布（分散は任意でよい。例えば1とすればよい）からサンプリングし、その採択/棄却を、上記のZ i,j(k)と同様にmin(1, exp J'/exp J)の確率にしたがって決める。
【００５７】
また、パラメータサンプリング部３２は、全てのパラメータについて更新すると、二次元配列Ｙ＾と、木構造サンプリング部３１により今回の反復で求められた木構造と、更新されたパラメータの各値とを用いて、上記（１３）式に従って、確率Ｊを算出し、更新されたパラメータの各値と共に確率Ｊを記憶部４に格納する。
【００５８】
終了判定部３３は、予め定められた終了条件を満たしたか否かを判定し、終了条件を満たしていない場合には、繰り返し回数を１つ増加させて、木構造サンプリング部３１、及びパラメータサンプリング部３２の各処理を繰り返す。本実施の形態では、生成するサンプルの個数の上限をあらかじめ指定しておき、当該サンプルの個数の上限を、繰り返し回数の上限として設定しておく。計算時間の許す限り出来るだけ大きな値を選ぶほど、サンプルは事後分布のより良い近似となるが、実用的には100個や数千個を選ぶのが一般的である。
【００５９】
終了判定部３３は、終了条件を満たしたと判定した場合には、記憶部４に記憶された、各反復で求められたパラメータC＾=(Cω,t,n)Ω×T×Nを全て出力部５により出力する。なお、パラメータC＾=(Cω,t,n)Ω×T×Nだけでなく、今までに生成した木構造Child＾={(Childleft(m))M,(Childright(m))M}、t＾={(tChildleft(m))M,(tChildright(m))M}のサンプル、各パラメータV＾=(Vn)N,τ＾=(τn)N,L＾=(Lk)K,Z＾=(Z(k)i,j)KK,K,β＾=(βk)Kの少なくとも１つを出力部５により出力してもよい。
【００６０】
上記のようなサンプリングを繰り返すことで、そこから生成されたサンプルは事後分布をよく近似したものとなる。
【００６１】
次に、本実施の形態に係る音響信号分析装置の作用について説明する。まず、分析対象の時系列信号として音響信号が音響信号分析装置に入力され、記憶部４に格納される。そして、音響信号分析装置において、図３に示す音響信号分析処理ルーチンが実行される。
【００６２】
まず、ステップ１００において、記憶部４から、各フレーム内の音響信号を読み込み、音響信号に対して、短時間フーリエ変換を用いた時間周波数分析を行った結果から、観測時間周波数成分Ｙω,tを各（ω，ｔ）の要素にもつ二次元配列Ｙ＾を生成して、記憶部４に記憶する。
【００６３】
そして、ステップ１０２において、乱数を用いて、パラメータC＾=(Cω,t,n)Ω×T×N, V＾=(Vn)N, t＾={(tChildleft(m))M,(tChildright(m))M}, τ＾=(τn)N,L＾=(Lk)K,Z＾=(Z(k)i,j)KK,K,β＾=(βk)Kの初期値を設定して、記憶部４に記憶する。
【００６４】
次にステップ１０４では、上記ステップ１０２で設定されたパラメータt＾={(tChildleft(m))M,(tChildright(m))M},L＾=(Lk)K,Z＾=(Z(k)i,j)KK,K、又は後述するステップ１０６で前回更新されたパラメータt＾={(tChildleft(m))M,(tChildright(m))M},L＾=(Lk)K,Z＾=(Z(k)i,j)KK,Kに基づいて、上記（１４）式に従って、四次元配列(pn,n',k,t)N×N×K×Tを再帰的に計算する。そして、計算された四次元配列(pn,n',k,t)N×N×K×Tを用いて、上記（１５）式に従って、木構造の各分割位置ｓの左右の子ノードのシンボル及び発音時刻をサンプリングにより決定し、Child＾={(Childleft(m))M,(Childright(m))M}、t＾={(tChildleft(m))M,(tChildright(m))M}を記憶部４に記憶する。
【００６５】
そして、ステップ１０６では、上記ステップ１０２で設定されたパラメータC＾=(Cω,t,n)Ω×T×N, V＾=(Vn)N, τ＾=(τn)N,L＾=(Lk)K,Z＾=(Z(k)i,j)KK,K,β＾=(βk)K、又は後述するステップ１０６で前回更新されたパラメータC＾=(Cω,t,n)Ω×T×N, V＾=(Vn)N, τ＾=(τn)N,L＾=(Lk)K,Z＾=(Z(k)i,j)KK,K,β＾=(βk)Kと、上記ステップ１００で生成された二次元行列Ｙ＾と、上記ステップ１０４で直前に得られたChild＾={(Childleft(m))M,(Childright(m))M}、t＾={(tChildleft(m))M,(tChildright(m))M}とに基づいて、サンプリングを行って、各パラメータの更新値C＾=(Cω,t,n)Ω×T×N, V＾=(Vn)N,, τ＾=(τn)N,L＾=(Lk)K,Z＾=(Z(k)i,j)KK,K,β＾=(βk)Kを決定し、記憶部４に記憶する。
【００６６】
次のステップ１０８では、所定の終了条件として、繰り返し回数が、上限に到達したか否かを判定し、繰り返し回数が上限に到達していない場合には、所定の終了条件が成立していないと判断して、上記ステップ１０４へ戻り、上記ステップ１０６で更新したパラメータを用いて、上記ステップ１０４〜ステップ１０６の処理を繰り返す。一方、繰り返し回数が上限に到達した場合には、所定の終了条件が成立したと判断し、ステップ１１０で、各反復の上記ステップ１０４、１０６で得られたパラメータC＾=(Cω,t,n)Ω×T×Nを出力部５により出力して、音響信号分析処理ルーチンを終了する。
【００６７】
次に、第２の実施の形態について説明する。なお、第１の実施の形態と同様の構成となる部分については、同一符号を付して説明を省略する。
【００６８】
第２の実施の形態では、各反復で得られたパラメータC＾=(Cω,t,n)Ω×T×Nのうち、何れかの反復で得られたパラメータC＾=(Cω,t,n)Ω×T×Nを選択して出力している点が、第１の実施の形態と異なっている。
【００６９】
図４に示すように、第２の実施の形態に係る音響信号分析装置の信号分析部３は、信号分析部３は、木構造サンプリング部３１と、パラメータサンプリング部３２と、終了判定部３３と、出力結果選択部３４とを備えている。
【００７０】
出力結果選択部３４は、記憶部４に記憶された、各反復で得られたパラメータC＾=(Cω,t,n)Ω×T×Nと、パラメータC＾=(Cω,t,n)Ω×T×Nを含む各パラメータの更新値に対応する確率Jとを入力とし、確率Jの値が最も高かった反復で得られたパラメータC＾=(Cω,t,n)Ω×T×Nを選択して、出力部５により出力する。
【００７１】
次に、音響音楽信号を分析対象として、上述した第１の実施の形態に係る手法を適用し、シミュレーション実験を行った結果について説明する。
【００７２】
ここで、入力信号として、非特許文献３（M. Goto, H. Hashiguchi, T. Nishimura, and R. Oka,"RWC music database: Popular, classical, and jazz music database" in ISMIR, 2002, pp. 287-288.）に記載されている音楽音響信号(RWC-MDB-C-2001 No.24A)をもともとの44.1kHzから16kHzに変換し、モノラルにしたものを使用した。振幅スペクトログラムを求める際には32msのハニング窓で16msごとにオーバーラップした短時間フーリエ変換を用いた。
【００７３】
図５（Ａ）に、楽器音の発音位置の正解としてのピアノロールを示し、図５（Ｂ）に、推定された各楽器音の音高および音量情報Σn VnOen,t-τnを示す。上記第１の実施の形態で説明した提案手法によれば、適切な数の楽器音と、それらの音高及び発音位置の情報を取得出来ていることが確認できた。
【００７４】
また、図６（Ａ）、（Ｂ）に、第１の実施の形態で説明した提案手法を適用した際にサンプリングにより生成された木構造の２つの例を示す。ｍ行は推定された木構造のｍ番目のノードのオンセット位置t（発音時刻）と長さL（シンボルの長さ）を表している。提案手法では、実際には2次元の木構造が抽出できているが、それらを描画するのは困難なため、上記図６（Ａ）、（Ｂ）においては各ノードの親子関係の情報を省略して表示している。例えば、上記図６（Ａ）の例において、28番目のノードは、同時出現する楽器音に分割する次元の分岐ルールによって36と37番目のノードに分割されている。上記図６（Ｂ）の例においては、32と33番目のノードは21番目のノードから時間軸方向に分割する次元の分岐ルールによって生成されている。また、上記図６（Ａ）、（Ｂ）において、各時刻のノードのうち、time frame軸に最も近いノードが、楽器音一音に対応している。
【００７５】
以上説明したように、本発明の実施の形態に係る音響信号分析装置によれば、音響信号の時系列データを頂点ノードとして、時間軸方向に音響信号を分割する次元と同時に発音する楽器音に分割する次元とにおいて、末端ノードが楽器音の一音に対応するように分割した木構造の各分割位置の子ノードのシンボル及び発音時刻のサンプリングと、Ｎ個の楽器音ｎの各々に対する時間周波数成分Ｃω,t,n、音量Ｖn、発音時刻に関するパラメータτn、各シンボルkに対する音の長さＬk、シンボルｋを持つ親ノードからシンボルｉを持つ子ノード及びシンボルｊを持つ子ノードに分岐させる重みを表わすパラメータＺi,j(k)、及び各シンボルｋの重みを表わすパラメータβkのサンプリングとを反復して行うことにより、楽器音一音単位に音響信号を分解すると共に、時間軸方向の分割と、楽器音が同時出現する場合とを考慮した構造解析を実現することができる。
【００７６】
また、複数の楽器音が混在する音響信号から、楽器音の識別を行うことなく各楽器音の一音ごとの信号を分離し、分離された一音一音が時間軸上のどの位置に配置されるかを精度良く分析することができる。
【００７７】
また、音響信号を入力として、その中の各音符情報とそれらの出現位置に関する構造を同時に推定する。その際に音符の出現位置に関して木構造の規則があることを仮定し、信号の構造解析を実現する。
【００７８】
なお、本発明は、上述した実施形態に限定されるものではなく、この発明の要旨を逸脱しない範囲内で様々な変形や応用が可能である。
【００７９】
例えば、上記の第２の実施の形態では、確率Ｊが最大となるときのパラメータＣ＾を一つ選択して、出力する場合を例に説明したが、これに限定されるものではない。分離結果の不確定さも含めて評価したい場合においては、反復毎に生成されたパラメータＣ＾のサンプルの中から無作為にいくつかのサンプルを取り出し、それらの頻度に従って分離結果の平均や分散、さらに高次のモーメントを評価するようにしてもよい。
【００８０】
また、上述の音響信号分析装置は、内部にコンピュータシステムを有しているが、「コンピュータシステム」は、ＷＷＷシステムを利用している場合であれば、ホームページ提供環境（あるいは表示環境）も含むものとする。
【００８１】
また、本願明細書中において、プログラムが予めインストールされている実施形態として説明したが、当該プログラムを、コンピュータ読み取り可能な記録媒体に格納して提供することも可能である。
【符号の説明】
【００８２】
１     時間周波数解析部
２     パラメータ初期値設定部
３     信号分析部
４     記憶部
５     出力部
３１   木構造サンプリング部
３２   パラメータサンプリング部
３３   終了判定部
３４   出力結果選択部
【図２】
2013195575.tif 000017
【図３】
2013195575.tif 000018
【図４】
2013195575.tif 000019
【図１】
2013195575.tif 000020
【図５】
2013195575.tif 000021
【図６】
2013195575.tif 000022