(57)【要約】【課題】音響信号の時系列データから、そこに含まれる音響イベントが、どのような音響要素で構成されているのかを精度よく推定することができるようにする。【解決手段】音響特徴量成分Ｃ＾を決定し、各時刻ｔに基底の状態が何れであるか、または新規の状態であることを表わすＺ^を決定し、基底の各状態ｋの基底スペクトルを表わすＨ^を決定し、予め与えられたＬ個の音響イベントｌの各々に対する一部の時刻ｔnの各々の音響イベントの有無をＸｌＮ^に基づいて、各時刻ｔの発音の有無を表すＵ^を決定し、遷移確率ａ（d）を決定すると共に遷移確率ｂ（d）を決定し、Ｌ個の音響イベントｌの各々に対するＤ個の音響要素ｄの各々の重みｗｌ＾を求め、Ｌ個の音響イベントｌの各々に対する精度パラメータαｌを算出し、予め定められた収束条件を満たすまで、各処理を繰り返し行い、Ｃ＾、Ｚ^、Ｈ^、Ｕ^及びｗｌ＾を出力する。【選択図】図４
2014164126.tif 000002
【特許請求の範囲】
【請求項１】
音響特徴量抽出手段、パラメータ初期値設定手段、音源分離手段、状態系列推論手段、基底スペクトル推論手段、アクティベーション推論手段、遷移確率推論手段、重みベクトル推論手段、精度パラメータ推論手段、及び収束判定手段を含む音響信号分析装置における音響信号分析方法であって、
前記音響特徴量抽出手段によって、複数の音響要素が混在する音響信号の時系列データを入力として、音響特徴量Ｙω,t（ωは周波数、ｔは時刻のインデックスである。）を要素にもつ二次元配列Ｙ＾を出力するステップと、
前記パラメータ初期値設定手段によって、予め定められたＤ個の音響要素ｄの各々における所定の基底の状態ｋの基底スペクトルを表わすパラメータＨω,d(k)を要素にもつ三次元配列Ｈ^、前記Ｄ個の音響要素の各々に対する各時刻ｔの発音の有無を表すパラメータＵd,t（＝｛０、１｝）を要素にもつ二次元配列Ｕ^、及び前記Ｄ個の音響要素ｄの各々に対して各時刻ｔに前記基底の状態が何れであるかを表わすパラメータＺｄ,tを要素にもつ二次元配列Ｚ^、予め定められたＬ個の音響イベントｌの各々に対する前記Ｄ個の音響要素ｄの各々の重みｗｌ，ｄを要素にもつ一次元配列ｗｌ＾、前記Ｌ個の音響イベントｌの各々に対する精度パラメータαｌの各々の初期値を設定すると共に、前記Ｄ個の音響要素ｄの各々について、前記パラメータＵd,tが０から１へ遷移する遷移確率ａd（ａ1＞ａ2＞・・・＞ａD-1＞ａD）、及び前記パラメータＵd,tが１から１へ遷移する遷移確率ｂdの各々の初期値を設定するステップと、
前記音源分離手段によって、（ω、ｔ、ｄ）の全ての組み合わせについて、前記二次元配列Ｚ^、前記三次元配列Ｈ^、及び前記二次元配列Ｕ^に基づいてパラメータｐω,t,dを算出し、算出した各パラメータｐω,t,dと、前記二次元配列Ｙ^とをパラメータとする多項分布に従って、Ｄ個の音響要素ｄの各々に対する音響特徴量成分Ｃω,t,ｄを要素にもつ三次元配列Ｃ＾をサンプリングにより決定するステップと、
前記状態系列推論手段によって、前記Ｄ個の音響要素ｄの各々について、前記二次元配列Ｚ^、前記三次元配列Ｃ^，及び前記二次元配列Ｕ^に基づいて、各時刻ｔに前記基底の状態が各状態ｋである事後確率、及び各時刻ｔに前記基底の状態が新規の状態である事後確率を算出し、算出した事後確率に従って、前記二次元配列Ｚ^をサンプリングにより決定するステップと、
前記基底スペクトル推論手段によって、（ｋ、ω、ｄ）の全ての組み合わせについて、前記二次元配列Ｚ^及び前記三次元配列Ｃ^に基づいて、パラメータφω,d(k)、ψω,d(k)を算出し、算出したパラメータφω,d(k)、ψω,d(k)をパラメータとする確率分布に従って、前記三次元配列Ｈ^をサンプリングにより決定するステップと、
前記アクティベーション推論手段によって、（ｄ、ｔ）の全ての組み合わせについて、予め与えられた前記Ｌ個の音響イベントｌの各々に対する一部の時刻ｔnの各々の音響イベントの有無を表すパラメータＸｌ,tｎ（＝｛０、１｝）を要素にもつ一次元配列ＸｌＮ^、前記三次元配列Ｃ^、前記二次元配列Ｚ^、前記三次元配列Ｈ^、及び前記遷移行列ａd，ｂd、及び前記一次元配列ｗｌ＾に基づいて、パラメータＵd,tの事後分布を算出し、算出したパラメータＵd,Tの事後分布に従って、Ｕd,Tをサンプリングにより決定し、ｔ＝Ｔ−１，・・・，１について、算出したパラメータＵd,ｔの事後分布と、前記遷移行列ａd，ｂdに基づく事後分布ｐ（Ｕd,ｔ＋１｜Ｕd,ｔ）との積に従って、Ｕd,ｔをサンプリングにより決定することにより、前記二次元配列Ｕ^を求めるステップと、
前記遷移確率推論手段によって、前記Ｄ個の音響要素ｄの各々について、前記遷移確率ａd-1、ａd＋1及び前記二次元配列Ｕ^に基づく遷移確率ａdの確率分布に従って、遷移確率ａdをサンプリングにより決定し、前記遷移確率ａDに基づく遷移確率ａdの確率分布に従って、遷移確率ａDをサンプリングにより決定すると共に、前記Ｄ個の音響要素ｄの各々について、前記二次元配列Ｕ^に基づく遷移確率ｂdの確率分布に従って、遷移確率ｂdをサンプリングにより決定するステップと、
前記重みベクトル推論手段によって、前記Ｌ個の音響イベントｌの各々について、前記一次元配列ＸｌＮ^、前記二次元配列Ｕ^、前記一次元配列ｗｌ＾、及び前記パラメータαｌに基づいて、事後確率が最大となるように前記一次元配列ｗｌ＾を求めるステップと、
前記精度パラメータ推論手段によって、前記Ｌ個の音響イベントｌの各々について、前記一次元配列ｗｌ＾に基づく前記パラメータαlの確率分布に従って、前記パラメータαｌをサンプリングにより決定するステップと、
前記収束判定手段によって、予め定められた収束条件を満たすまで、前記音源分離手段、前記状態系列推論手段、前記基底スペクトル推論手段、前記アクティベーション推論手段、前記遷移確率推論手段、前記重みベクトル推論手段、及び前記精度パラメータ推論手段による各処理を繰り返し行い、前記三次元配列Ｃ＾、前記二次元配列Ｚ^、前記三次元配列Ｈ^、前記二次元配列Ｕ^、及び前記Ｌ個の音響イベントｌの各々の前記一次元配列ｗｌ＾を出力するステップと、
を含む音響信号分析方法。
【請求項２】
音響イベントラベル推論手段によって、前記Ｌ個の音響イベントｌの各々について、前記二次元配列Ｕ^、前記一次元配列ＸｌＮ^、前記精度パラメータαｌ、及び前記一次元配列ｗｌ＾に基づいて、各時刻ｔのパラメータＸｌ,tの予測確率を算出し、前記予測確率に基づいて、各時刻ｔのパラメータＸｌ,tの値を求めるステップを更に含む請求項１記載の音響信号分析方法。
【請求項３】
前記初期値設定手段によって、初期値を設定するステップは、
前記予め定められたＤ個の音響要素ｄの各々における基底の状態数を１とし、全てのパラメータＨω,d(1)が非負値であり、全てのパラメータＵd,tが非負値である、という条件の下で、前記パラメータＨω,d(1)と前記パラメータＵd,tとの積を、全てのｄについて足し合わせたモデルについて、前記音響特徴量抽出手段によって出力された音響特徴量Ｙω,tと前記モデルとの距離を表わした目的関数の値が小さくなるように、前記パラメータＨω,d(1)及び前記パラメータＵd,tの各々を推定し、前記推定された前記パラメータＵd,tの各々を二値化することにより、前記三次元配列Ｈ^及び前記二次元配列Ｕ^の初期値を設定すると共に、前記パラメータＺｄ,tの各々を１とした二次元配列Ｚ^の初期値を設定する請求項１又は２記載の音響信号分析方法。
【請求項４】
複数の音響要素が混在する音響信号の時系列データを入力として、音響特徴量Ｙω,t（ωは周波数、ｔは時刻のインデックスである。）を要素にもつ二次元配列Ｙ＾を出力する音響特徴量抽出手段と、
予め定められたＤ個の音響要素ｄの各々における所定の基底の状態ｋの基底スペクトルを表わすパラメータＨω,d(k)を要素にもつ三次元配列Ｈ^、前記Ｄ個の音響要素の各々に対する各時刻ｔの発音の有無を表すパラメータＵd,t（＝｛０、１｝）を要素にもつ二次元配列Ｕ^、及び前記Ｄ個の音響要素ｄの各々に対して各時刻ｔに前記基底の状態が何れであるかを表わすパラメータＺｄ,tを要素にもつ二次元配列Ｚ^、予め定められたＬ個の音響イベントｌの各々に対する前記Ｄ個の音響要素ｄの各々の重みｗｌ，ｄを要素にもつ一次元配列ｗｌ＾、前記Ｌ個の音響イベントｌの各々に対する精度パラメータαｌの各々の初期値を設定すると共に、前記Ｄ個の音響要素ｄの各々について、前記パラメータＵd,tが０から１へ遷移する遷移確率ａd（ａ1＞ａ2＞・・・＞ａD-1＞ａD）、及び前記パラメータＵd,tが１から１へ遷移する遷移確率ｂdの各々の初期値を設定するパラメータ初期値設定手段と、
（ω、ｔ、ｄ）の全ての組み合わせについて、前記二次元配列Ｚ^、前記三次元配列Ｈ^、及び前記二次元配列Ｕ^に基づいてパラメータｐω,t,dを算出し、算出した各パラメータｐω,t,dと、前記二次元配列Ｙ^とをパラメータとする多項分布に従って、Ｄ個の音響要素ｄの各々に対する音響特徴量成分Ｃω,t,ｄを要素にもつ三次元配列Ｃ＾をサンプリングにより決定する音源分離手段と、
前記Ｄ個の音響要素ｄの各々について、前記二次元配列Ｚ^、前記三次元配列Ｃ^，及び前記二次元配列Ｕ^に基づいて、各時刻ｔに前記基底の状態が各状態ｋである事後確率、及び各時刻ｔに前記基底の状態が新規の状態である事後確率を算出し、算出した事後確率に従って、前記二次元配列Ｚ^をサンプリングにより決定する状態系列推論手段と、
（ｋ、ω、ｄ）の全ての組み合わせについて、前記二次元配列Ｚ^及び前記三次元配列Ｃ^に基づいて、パラメータφω,d(k)、ψω,d(k)を算出し、算出したパラメータφω,d(k)、ψω,d(k)をパラメータとする確率分布に従って、前記三次元配列Ｈ^をサンプリングにより決定する基底スペクトル推論手段と、
（ｄ、ｔ）の全ての組み合わせについて、予め与えられた前記Ｌ個の音響イベントｌの各々に対する一部の時刻ｔnの各々の音響イベントの有無を表すパラメータＸｌ,tｎ（＝｛０、１｝）を要素にもつ一次元配列ＸｌＮ^、前記三次元配列Ｃ^、前記二次元配列Ｚ^、前記三次元配列Ｈ^、及び前記遷移行列ａd，ｂd、及び前記一次元配列ｗｌ＾に基づいて、パラメータＵd,tの事後分布を算出し、算出したパラメータＵd,Tの事後分布に従って、Ｕd,Tをサンプリングにより決定し、ｔ＝Ｔ−１，・・・，１について、算出したパラメータＵd,ｔの事後分布と、前記遷移行列ａd，ｂdに基づく事後分布ｐ（Ｕd,ｔ＋１｜Ｕd,ｔ）との積に従って、Ｕd,ｔをサンプリングにより決定することにより、前記二次元配列Ｕ^を求めるアクティベーション推論手段と、
前記Ｄ個の音響要素ｄの各々について、前記遷移確率ａd-1、ａd＋1及び前記二次元配列Ｕ^に基づく遷移確率ａdの確率分布に従って、遷移確率ａdをサンプリングにより決定し、前記遷移確率ａDに基づく遷移確率ａdの確率分布に従って、遷移確率ａDをサンプリングにより決定すると共に、前記Ｄ個の音響要素ｄの各々について、前記二次元配列Ｕ^に基づく遷移確率ｂdの確率分布に従って、遷移確率ｂdをサンプリングにより決定する遷移確率推論手段と、
前記Ｌ個の音響イベントｌの各々について、前記一次元配列ＸｌＮ^、前記二次元配列Ｕ^、前記一次元配列ｗｌ＾、及び前記パラメータαｌに基づいて、事後確率が最大となるように前記一次元配列ｗｌ＾を求める重みベクトル推論手段と、
前記Ｌ個の音響イベントｌの各々について、前記一次元配列ｗｌ＾に基づく前記パラメータαlの確率分布に従って、前記パラメータαｌをサンプリングにより決定する精度パラメータ推論手段と、
予め定められた収束条件を満たすまで、前記音源分離手段、前記状態系列推論手段、前記基底スペクトル推論手段、前記アクティベーション推論手段、前記遷移確率推論手段、前記重みベクトル推論手段、及び前記精度パラメータ推論手段による各処理を繰り返し行い、前記三次元配列Ｃ＾、前記二次元配列Ｚ^、前記三次元配列Ｈ^、前記二次元配列Ｕ^、及び前記Ｌ個の音響イベントｌの各々の前記一次元配列ｗｌ＾を出力する収束判定手段と、
を含む音響信号分析装置。
【請求項５】
請求項１〜請求項３の何れか１項記載の音響信号分析方法を構成する各ステップをコンピュータに実行させるためのプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、音響信号分析方法、装置、及びプログラムに係り、特に、音響信号の時系列から、信号パラメータを分析する音響信号分析方法、装置、及びプログラムに関する。
【背景技術】
【０００２】
膨大な音や映像のメディアデータが身の回りにあふれる中、これらのデータを自在に検索して活用するためには、付随するテキストデータに頼るだけではなく、それぞれの中身を表す情報を、音や映像自体から自動的に引き出す技術が必要不可欠である。
【０００３】
従来は、教師あり学習の下で、音響信号に含まれる音響イベントを自動的にラベリングする技術が一般的であった。すなわち、学習データとして、ラベル付けされた音響イベントの音源信号を事前に用意し、それを周波数分析して得られる音響的特徴と統計モデルを用いて、未知の音響信号の時系列データにラベル付けを行った（例えば、非特許文献１、非特許文献２、非特許文献３）。統計モデルでは、例えば、ガウス混合モデル（ＧＭＭ） を用いる場合、各音響イベントを表現しうる音響的特徴の頻度分布が学習される。隠れマルコフモデル（ＨＭＭ）を用いる場合、音響的特徴の統計的な時間遷移が学習され、これを用いて未知の音響信号の時系列データにラベル付けを行う。これらのラベル付けの多段処理技術（音響信号が無音なのかそうでないか、音楽か話声かを段階的に識別する）も提案されている（例えば、非特許文献４）。また、音の重ね合わせを考慮して、音響信号をあらかじめ、非負値行列分解（ＮＭＦ）によって音源分離し、分離信号を学習データに用いる技術も提案されている（例えば、非特許文献５）。
【先行技術文献】
【非特許文献】
【０００４】
【非特許文献１】J. Saunders、“ Real-time discrimination of broadcast speech/music, ”、in Proc. ICASSP、1996．
【非特許文献２】T. Butko and C. Nadeu、“ Audio segmentation of broadcast news in the Albayzin-2010 evaluation: overview, results, and discussion, ”、EURASIP Journal on Audio, Speech, and Music Processing、2011．
【非特許文献３】A. Mesaros, T. Heittola, A. Eronen and T. Virtanen、“ Acoustic event detection in real life recordings,”、in Proc. EUSIPCO、2010．
【非特許文献４】T. Butko and C. Nadeu、“ Audio segmentation of broadcast news: A hierarchical system with feature selection for the Albayzin-2010 evaluation, ”、in Proc.2 ICASSP、2011．
【非特許文献５】T. Heittola, A. Mesaros, T. Virtanen and A. Eronen、“ Sound Event Detection in Multisource Environments Using Source Separation, ”、in Proc. CHiME、2011．
【発明の概要】
【発明が解決しようとする課題】
【０００５】
しかし従来の技術では、以下の２つの課題があった。１つ目の課題は音響イベントの重なりがこれまで十分に考慮されなかった点である。多くの研究では、イベントの重ね合わせは無視し、多重音の音響信号から特徴を抽出して音響モデルを学習して、各フレームで最も顕著なイベントだけを出力した。音響イベントがスパースに現れる音環境では十分な性能が得られるものの、音響イベントが豊富に含まれる環境を対象とすると検出が難しい。また、あらかじめ音源分離を行い分離信号を用いれば、音響イベントの重なりを考慮できるが、音環境に合わせて分離信号のトラック数を手動で調整する必要があった。
【０００６】
２つ目の課題は、数千時間の書き起こしデータを用いる音声認識に比べ、音響イベントの音響的特徴を学習するためのラベル付データベースが少ない点である。人手でラベル付された音響信号データを構築するためにはコストがかかる。また、そもそも観測される音響信号に、どのような音響イベントが含まれているか事前知識がないことも多く、音響イベントの総数や音響イベントを表現するための音響的特徴はモデルが大量の音響信号データから自動的に決定してくれることが理想的である。
【０００７】
本発明は、上記の事情を考慮してなされたもので、音響信号の時系列データから、そこに含まれる音響イベントが、どのような音響要素で構成されているのかを精度よく推定することができる音響信号分析方法、装置、及びプログラムを提供することを目的とする。
【課題を解決するための手段】
【０００８】
上記の目的を達成するために本発明に係る音響信号分析方法は、音響特徴量抽出手段、パラメータ初期値設定手段、音源分離手段、状態系列推論手段、基底スペクトル推論手段、アクティベーション推論手段、遷移確率推論手段、重みベクトル推論手段、精度パラメータ推論手段、及び収束判定手段を含む音響信号分析装置における音響信号分析方法であって、前記音響特徴量抽出手段によって、複数の音響要素が混在する音響信号の時系列データを入力として、音響特徴量Ｙω,t（ωは周波数、ｔは時刻のインデックスである。）を要素にもつ二次元配列Ｙ＾を出力するステップと、前記パラメータ初期値設定手段によって、予め定められたＤ個の音響要素ｄの各々における所定の基底の状態ｋの基底スペクトルを表わすパラメータＨω,d(k)を要素にもつ三次元配列Ｈ^、前記Ｄ個の音響要素の各々に対する各時刻ｔの発音の有無を表すパラメータＵd,t（＝｛０、１｝）を要素にもつ二次元配列Ｕ^、及び前記Ｄ個の音響要素ｄの各々に対して各時刻ｔに前記基底の状態が何れであるかを表わすパラメータＺｄ,tを要素にもつ二次元配列Ｚ^、予め定められたＬ個の音響イベントｌの各々に対する前記Ｄ個の音響要素ｄの各々の重みｗｌ，ｄを要素にもつ一次元配列ｗｌ＾、前記Ｌ個の音響イベントｌの各々に対する精度パラメータαｌの各々の初期値を設定すると共に、前記Ｄ個の音響要素ｄの各々について、前記パラメータＵd,tが０から１へ遷移する遷移確率ａd（ａ1＞ａ2＞・・・＞ａD-1＞ａD）、及び前記パラメータＵd,tが１から１へ遷移する遷移確率ｂdの各々の初期値を設定するステップと、前記音源分離手段によって、（ω、ｔ、ｄ）の全ての組み合わせについて、前記二次元配列Ｚ^、前記三次元配列Ｈ^、及び前記二次元配列Ｕ^に基づいてパラメータｐω,t,dを算出し、算出した各パラメータｐω,t,dと、前記二次元配列Ｙ^とをパラメータとする多項分布に従って、Ｄ個の音響要素ｄの各々に対する音響特徴量成分Ｃω,t,ｄを要素にもつ三次元配列Ｃ＾をサンプリングにより決定するステップと、前記状態系列推論手段によって、前記Ｄ個の音響要素ｄの各々について、前記二次元配列Ｚ^、前記三次元配列Ｃ^，及び前記二次元配列Ｕ^に基づいて、各時刻ｔに前記基底の状態が各状態ｋである事後確率、及び各時刻ｔに前記基底の状態が新規の状態である事後確率を算出し、算出した事後確率に従って、前記二次元配列Ｚ^をサンプリングにより決定するステップと、前記基底スペクトル推論手段によって、（ｋ、ω、ｄ）の全ての組み合わせについて、前記二次元配列Ｚ^及び前記三次元配列Ｃ^に基づいて、パラメータφω,d(k)、ψω,d(k)を算出し、算出したパラメータφω,d(k)、ψω,d(k)をパラメータとする確率分布に従って、前記三次元配列Ｈ^をサンプリングにより決定するステップと、前記アクティベーション推論手段によって、（ｄ、ｔ）の全ての組み合わせについて、予め与えられた前記Ｌ個の音響イベントｌの各々に対する一部の時刻ｔnの各々の音響イベントの有無を表すパラメータＸｌ,tｎ（＝｛０、１｝）を要素にもつ一次元配列ＸｌＮ^、前記三次元配列Ｃ^、前記二次元配列Ｚ^、前記三次元配列Ｈ^、及び前記遷移行列ａd，ｂd、及び前記一次元配列ｗｌ＾に基づいて、パラメータＵd,tの事後分布を算出し、算出したパラメータＵd,Tの事後分布に従って、Ｕd,Tをサンプリングにより決定し、ｔ＝Ｔ−１，・・・，１について、算出したパラメータＵd,ｔの事後分布と、前記遷移行列ａd，ｂdに基づく事後分布ｐ（Ｕd,ｔ＋１｜Ｕd,ｔ）との積に従って、Ｕd,ｔをサンプリングにより決定することにより、前記二次元配列Ｕ^を求めるステップと、前記遷移確率推論手段によって、前記Ｄ個の音響要素ｄの各々について、前記遷移確率ａd-1、ａd＋1及び前記二次元配列Ｕ^に基づく遷移確率ａdの確率分布に従って、遷移確率ａdをサンプリングにより決定し、前記遷移確率ａDに基づく遷移確率ａdの確率分布に従って、遷移確率ａDをサンプリングにより決定すると共に、前記Ｄ個の音響要素ｄの各々について、前記二次元配列Ｕ^に基づく遷移確率ｂdの確率分布に従って、遷移確率ｂdをサンプリングにより決定するステップと、前記重みベクトル推論手段によって、前記Ｌ個の音響イベントｌの各々について、前記一次元配列ＸｌＮ^、前記二次元配列Ｕ^、前記一次元配列ｗｌ＾、及び前記パラメータαｌに基づいて、事後確率が最大となるように前記一次元配列ｗｌ＾を求めるステップと、前記精度パラメータ推論手段によって、前記Ｌ個の音響イベントｌの各々について、前記一次元配列ｗｌ＾に基づく前記パラメータαlの確率分布に従って、前記パラメータαｌをサンプリングにより決定するステップと、前記収束判定手段によって、予め定められた収束条件を満たすまで、前記音源分離手段、前記状態系列推論手段、前記基底スペクトル推論手段、前記アクティベーション推論手段、前記遷移確率推論手段、前記重みベクトル推論手段、及び前記精度パラメータ推論手段による各処理を繰り返し行い、前記三次元配列Ｃ＾、前記二次元配列Ｚ^、前記三次元配列Ｈ^、前記二次元配列Ｕ^、及び前記Ｌ個の音響イベントｌの各々の前記一次元配列ｗｌ＾を出力するステップと、を含む。
【０００９】
本発明に係る音響信号分析装置は、  複数の音響要素が混在する音響信号の時系列データを入力として、音響特徴量Ｙω,t（ωは周波数、ｔは時刻のインデックスである。）を要素にもつ二次元配列Ｙ＾を出力する音響特徴量抽出手段と、予め定められたＤ個の音響要素ｄの各々における所定の基底の状態ｋの基底スペクトルを表わすパラメータＨω,d(k)を要素にもつ三次元配列Ｈ^、前記Ｄ個の音響要素の各々に対する各時刻ｔの発音の有無を表すパラメータＵd,t（＝｛０、１｝）を要素にもつ二次元配列Ｕ^、及び前記Ｄ個の音響要素ｄの各々に対して各時刻ｔに前記基底の状態が何れであるかを表わすパラメータＺｄ,tを要素にもつ二次元配列Ｚ^、予め定められたＬ個の音響イベントｌの各々に対する前記Ｄ個の音響要素ｄの各々の重みｗｌ，ｄを要素にもつ一次元配列ｗｌ＾、前記Ｌ個の音響イベントｌの各々に対する精度パラメータαｌの各々の初期値を設定すると共に、前記Ｄ個の音響要素ｄの各々について、前記パラメータＵd,tが０から１へ遷移する遷移確率ａd（ａ1＞ａ2＞・・・＞ａD-1＞ａD）、及び前記パラメータＵd,tが１から１へ遷移する遷移確率ｂdの各々の初期値を設定するパラメータ初期値設定手段と、（ω、ｔ、ｄ）の全ての組み合わせについて、前記二次元配列Ｚ^、前記三次元配列Ｈ^、及び前記二次元配列Ｕ^に基づいてパラメータｐω,t,dを算出し、算出した各パラメータｐω,t,dと、前記二次元配列Ｙ^とをパラメータとする多項分布に従って、Ｄ個の音響要素ｄの各々に対する音響特徴量成分Ｃω,t,ｄを要素にもつ三次元配列Ｃ＾をサンプリングにより決定する音源分離手段と、前記Ｄ個の音響要素ｄの各々について、前記二次元配列Ｚ^、前記三次元配列Ｃ^，及び前記二次元配列Ｕ^に基づいて、各時刻ｔに前記基底の状態が各状態ｋである事後確率、及び各時刻ｔに前記基底の状態が新規の状態である事後確率を算出し、算出した事後確率に従って、前記二次元配列Ｚ^をサンプリングにより決定する状態系列推論手段と、（ｋ、ω、ｄ）の全ての組み合わせについて、前記二次元配列Ｚ^及び前記三次元配列Ｃ^に基づいて、パラメータφω,d(k)、ψω,d(k)を算出し、算出したパラメータφω,d(k)、ψω,d(k)をパラメータとする確率分布に従って、前記三次元配列Ｈ^をサンプリングにより決定する基底スペクトル推論手段と、（ｄ、ｔ）の全ての組み合わせについて、予め与えられた前記Ｌ個の音響イベントｌの各々に対する一部の時刻ｔnの各々の音響イベントの有無を表すパラメータＸｌ,tｎ（＝｛０、１｝）を要素にもつ一次元配列ＸｌＮ^、前記三次元配列Ｃ^、前記二次元配列Ｚ^、前記三次元配列Ｈ^、及び前記遷移行列ａd，ｂd、及び前記一次元配列ｗｌ＾に基づいて、パラメータＵd,tの事後分布を算出し、算出したパラメータＵd,Tの事後分布に従って、Ｕd,Tをサンプリングにより決定し、ｔ＝Ｔ−１，・・・，１について、算出したパラメータＵd,ｔの事後分布と、前記遷移行列ａd，ｂdに基づく事後分布ｐ（Ｕd,ｔ＋１｜Ｕd,ｔ）との積に従って、Ｕd,ｔをサンプリングにより決定することにより、前記二次元配列Ｕ^を求めるアクティベーション推論手段と、前記Ｄ個の音響要素ｄの各々について、前記遷移確率ａd-1、ａd＋1及び前記二次元配列Ｕ^に基づく遷移確率ａdの確率分布に従って、遷移確率ａdをサンプリングにより決定し、前記遷移確率ａDに基づく遷移確率ａdの確率分布に従って、遷移確率ａDをサンプリングにより決定すると共に、前記Ｄ個の音響要素ｄの各々について、前記二次元配列Ｕ^に基づく遷移確率ｂdの確率分布に従って、遷移確率ｂdをサンプリングにより決定する遷移確率推論手段と、前記Ｌ個の音響イベントｌの各々について、前記一次元配列ＸｌＮ^、前記二次元配列Ｕ^、前記一次元配列ｗｌ＾、及び前記パラメータαｌに基づいて、事後確率が最大となるように前記一次元配列ｗｌ＾を求める重みベクトル推論手段と、前記Ｌ個の音響イベントｌの各々について、前記一次元配列ｗｌ＾に基づく前記パラメータαlの確率分布に従って、前記パラメータαｌをサンプリングにより決定する精度パラメータ推論手段と、予め定められた収束条件を満たすまで、前記音源分離手段、前記状態系列推論手段、前記基底スペクトル推論手段、前記アクティベーション推論手段、前記遷移確率推論手段、前記重みベクトル推論手段、及び前記精度パラメータ推論手段による各処理を繰り返し行い、前記三次元配列Ｃ＾、前記二次元配列Ｚ^、前記三次元配列Ｈ^、前記二次元配列Ｕ^、及び前記Ｌ個の音響イベントｌの各々の前記一次元配列ｗｌ＾を出力する収束判定手段と、を含んで構成されている。
【００１０】
本発明に係るプログラムは、上記の音響信号分析方法を構成する各ステップをコンピュータに実行させるためのプログラムである。
【発明の効果】
【００１１】
以上説明したように、本発明の音響信号分析方法、装置、及びプログラムによれば、音響要素ｄの各々に対する音響特徴量成分Ｃω,t,ｄを要素にもつＣ＾をサンプリングにより決定し、音響要素ｄの各々に対して各時刻ｔに基底の状態が何れであるか、または新規の状態であることを表わすパラメータＺｄ,tを要素にもつＺ^をサンプリングにより決定し、音響要素ｄの各々における基底の各状態ｋの基底スペクトルを表わすパラメータＨω,d(k)を要素にもつＨ^をサンプリングにより決定し、予め与えられたＬ個の音響イベントｌの各々に対する一部の時刻ｔnの各々の音響イベントの有無を表すパラメータＸｌ,tｎ（＝｛０、１｝）を要素にもつＸｌＮ^を用いて、音響要素ｄの各々に対する各時刻ｔの発音の有無を表すパラメータＵd,tを要素にもつＵ^をサンプリングにより決定し、音響要素ｄの各々について、パラメータＵd,tが０から１へ遷移する遷移確率ａ（d）をサンプリングにより決定すると共にパラメータＵd,tが１から１へ遷移する遷移確率ｂ（d）をサンプリングにより決定し、予め与えられたＸｌＮ^を用いて、音響イベントｌの各々に対するＤ個の音響要素ｄの各々の重みｗｌ，ｄを要素にもつｗｌ＾を求め、音響イベントｌの各々に対する精度パラメータαｌを算出することを、予め定められた収束条件を満たすまで繰り返し行い、Ｃ＾、Ｚ^、Ｈ^、Ｕ^及びｗｌ＾を出力することにより、音響信号の時系列データから、そこに含まれる音響イベントが、どのような音響要素で構成されているのかを精度よく推定することができる、という効果が得られる。
【図面の簡単な説明】
【００１２】
【図１】音響信号への音響イベントのマルチラベリングを示すイメージ図である。
【図２】音響イベントの音響要素を構成する基底スペクトル集合、状態系列集合、及びアクティベーション集合を示すイメージ図である。
【図３】本発明の実施の形態に係る提案モデルのグラフィカル表現図である。
【図４】本発明の実施の形態に係る音響信号分析装置の構成を示す概略図である。
【図５】本発明の実施の形態に係る音響信号分析装置における音響信号分析処理ルーチンの内容を示すフローチャートである。
【図６】本発明の実施の形態の実験結果を示す図である。
【図７】本発明の実施の形態における学習データの長さに対する実験結果を示す図である。
【図８】本発明の実施の形態と従来技術との性能を比較した実験結果を示す図である。
【図９】アクティベーション行列に推定される未知の音響イベントを示す図である。
【発明を実施するための形態】
【００１３】
以下、図面を参照して本発明の実施の形態を詳細に説明する。以下では、まず本発明の実施の形態の概要及び原理について説明する。
【００１４】
＜発明の概要＞
本発明の実施の形態は、テレビやラジオ放送、ポッドキャスト、動画投稿サイトにおける音響信号の時系列データに対して、そこに含まれる、人に認識されうる音の事象“音響イベント”（話声や歌声、笑い声やあいづちをはじめ、動物の鳴き声、楽器音、環境音、効果音などの音のカテゴリを指す）を書き起こすアルゴリズムに関する。具体的には、音響信号を周波数分析して得られるスペクトログラムをフィルタバンク処理し、その出力値行列に含まれる、様々な音響要素（音声における音素、音楽における楽器音など）の音響的特徴を表現するための基底スペクトル集合、各音響要素の基底スペクトルの状態遷移を表現する状態系列集合、そして各音響要素の発音区間（鳴っているか否か、ＯＮ／ＯＦＦの状態からなる）を表現するアクティベーション集合を抽出し、抽出されたアクティベーションのベイズロジスティック回帰に基づいて、フレームごとに音響イベントをラベル付けする信号解析装置、信号解析方法及び信号解析プログラムに関する。図１は、音響信号から音響イベントを書き起こすための概略図を示す。また、図２は、フィルタバンク出力値を分解して得られる音響要素の基底スペクトル集合、状態系列集合、アクティベーション集合の概略図を示す。
【００１５】
非特許文献６（大石康智、持橋大地、松井知子、中野允裕、亀岡弘和、泉谷知範、柏野邦夫、「ノンパラメトリックベイズアプローチに基づく音響イベント検出」、電子情報通信学会技術研究報告 Vol.112,No.111、2012年6月22日、p.37-42）では、機械学習分野で注目されるノンパラメトリックベイズ法を利用して、大規模な音響信号の時系列データから、そこに含まれる全ての音響要素の音響的特徴とその発音区間を自律的に学習させるフレームワークを提供した（音響イベントはこの音響要素を組み合わせることによって構成されると想定する）。具体的には、音響信号に含まれる音響要素の数（上記図２に示すＤに相当する）が無限個の可能性を持つと仮定し、各音響要素の継続時間を表現するためにマルコフ性を導入したMarkov Indian Buffet Process を用いて、各音響要素の音響的特徴およびそれらのアクティベーションを確率的に推論する。本実施の形態では、推定された音響要素のアクティベーションを利用して、ベイズロジスティック回帰によって、各フレームに音響イベントをラベル付けすることを提供する。
【００１６】
＜原理＞
[１.  音響イベント転写モデル]
音響イベントの重ねあわせをモデル化するために、可変基底型ＮＭＦ（例えば、参考文献１（A. Ozerov, C. F´evotte, and M. Charbit,“ Factorial scaled hidden Markov model for polyphonic audio representation and source separation, ”in Proc. WASPAA 2009. ）や、参考文献２（M. Nakano, J. Le Roux, H. Kameoka, N. Ono and S. Sagayama,“ Infinite-state spectrum model for music signal analysis, ”in Proc. ICASSP 2011.）を参照）を利用する。これは各音響要素の基底スペクトルが時間にともなって変化するように拡張したＮＭＦである。各音響要素の音色を特徴付けるため、振幅スペクトログラムをメルフィルタバンク処理して得られた出力値行列である音響特徴量行列Ｙ＾＝（Ｙω，ｔ）Ω×T∈Ｒ≧0,Ω×Ｔに、このＮＭＦを適用する。なお、記号に付された「＾」は、当該記号が行列または多次元配列またはベクトルである。ここで、ω＝１,...,Ωはメルフィルタバンクの中心周波数のインデックス、ｔ＝１,...,Ｔは時間のインデックス、ｄ＝１,...,Ｄは音響要素のインデックスを表す。このとき、各音響要素の基底スペクトルは時刻t に、ある一つの状態Ｚｄ,ｔ∈Ｎを取ると見なし、
【００１７】
【数１】
2014164126.tif 000003
【００１８】
と表現する。これは図２に示すような要素からなる生成モデルである。まず、Ｈｄ＾＝（Ｈ（ｋ）ω，ｄ）Ω×Ｋｄは音響要素ｄのＫｄ個の基底スペクトルを意味する。Ｕ＾＝（Ｕｄ，ｔ）Ｄ×Tは、音響要素のｏｎ／ｏｆｆを表現する２値系列からなるアクティベーション集合を意味する。このとき、音響信号の音量に関する情報はすべて、基底スペクトルに含まれて表現されると想定する。そして、Ｃｄ＾＝（Ｃω，t，d）Ω×Ｔは、音響要素ｄのメルフィルタバンク出力値成分である音響特徴量成分に相当し、Ｃω，t，dはＨ（ｋ）ω，ｄ，Ｕｄ，ｔ，Ｚｄ,ｔをパラメータとするポアソン分布から生成される。最終的に音響信号はこれらの音響要素の和で表現される。先に述べたように、音響イベントはこれらの音響要素の組み合わせで表現されることを想定している。
さらに、Ｕ＾とＺ＾の事前分布としてｍＩＢＰとＣＲＰを導入し、無限個の音響要素、無限個の基底スペクトルを仮定するモデルへと拡張する。また、半教師あり学習の下で、ベイズロジスティック回帰を利用して、音響要素の組み合わせに音響イベントをラベル付けるアルゴリズムを導入する。以下にその導入方法を説明する。
【００１９】
[１．１  ｍＩＢＰによって生成されるアクティベーション行列]
Markov Indian buffet process(ｍＩＢＰ)と呼ばれるノンパラメトリックベイジアン因子モデル（参考文献３（J. V. Gael et al., “ The infinite factorial hidden Markov model, ”in Proc. NIPS 2008）を参照）は、時間的な依存性を考慮するようにIndian buffet process (IBP)（参考文献４（T. L. Griffiths et al.,“ Infinite latent feature models and the Indian buffet process, ”in Proc. NIPS 2006.）を参照）を拡張したものである。これは各フレームにおける音響要素がｏｎかｏｆｆかどうかをモデル化するためにアクティベーション行列Ｕ＾に対する事前分布を定義する。次の２つの特性を満たす。
１.行数（音響要素の数に相当する）は任意に大きいサイズを想定する
２.各列（フレームに相当する）のｏｎ／ｏｆｆは、音響要素の継続時間（発音時間）を表現するよう、マルコフ過程に従って生起されるこの構成方法が時系列データに対する因子表現を学習することを可能にする。このｍＩＢＰをアクティベーション行列Ｕ＾の事前分布に導入する。具体的には、音響要素の発音時間を表現するために、各音響要素ごとに遷移行列
【００２０】
【数２】
2014164126.tif 000004
【００２１】
を用意する。つまり、０→０の遷移確率は１−ａｄ，０→１の遷移確率はａｄ，１→０の遷移確率は１−ｂｄ，１→１の遷移確率はｂｄとし、これらの状態遷移によって、０と１からなる音響要素ｄのアクティベーション系列を生成する。Ｕｄ，０＝０とすると、
【００２２】
【数３】
2014164126.tif 000005
【００２３】
のようなベルヌーイ分布から生成される確率変数として記述できる。このとき、すべての音響要素のアクティベーションからなるアクティベーション行列Ｕ＾は
【００２４】
【数４】
2014164126.tif 000006
【００２５】
に従う。ここで、ａ＾＝{ａ１，．．．，ａＤ}、ｂ＾＝{ｂ１，．．．，ｂＤ}とし、ｃｄ００、ｃｄ０１、ｃｄ１０、ｃｄ１１はそれぞれ０→０、０→１、１→０、１→１に遷移する回数とする。さらに、ａｄとｂｄの事前分布をそれらの共役性から、
【００２６】
【数５】
2014164126.tif 000007
【００２７】
とする。式（４）における有限モデルをＤ→∞として無限モデルに拡張するために、stick breaking construction（参考文献５（Y. W. Teh et al.,“ Stick-breaking construction for the Indian buffet process, ” in Proc. NIPS 2007）参照）を適用する。ここで、ａ＾をａ(1)＞ａ(2)＞．．．＞ａ(D)のように順序付けると、Ｄ→∞におけるａ（ｄ）の生成過程は以下の法則に従う。
【００２８】
【数６】
2014164126.tif 000008
【００２９】
一方、ｂｄはベータ分布Ｂｅｔａ（θｂ（０），θｂ（１））から独立にサンプルされる。したがって、ｄ番目に大きいａｄに対応する変数をｂ（ｄ）とすると、ｂ（ｄ）〜Ｂｅｔａ（θｂ（０），θｂ（１））となる。ｍＩＢＰでは、頻繁にアクティベートされる音響要素ほど、小さいインデックスｄに割り当てられる。
【００３０】
[１．２  ＣＲＰによって生成される時変なスペクトル特性]
音響イベントを表現するために必要な、基底スペクトルの状態数は音響イベントごとに異なる。例えば、ピアノの単音であれば、“attack”、“decay”、“sustain”、“release”と呼ばれる状態があるため、４つの基底スペクトルを用いて表現されるだろう。同様のことが音声を構成する音素にも言える。ドアの開閉音や食器音は、楽器の単音に比べて突発音であるため、基底スペクトルの状態数は少ないかもしれない。交通騒音（ノイズ）は常に1 つの基底スペクトルで表現されるかもしれない。このように、基底スペクトルの状態数は固定するのではなく、音響信号から自動的に決定されることが望ましい。本発明者らは可変基底型ＮＭＦにディリクレ過程を導入した無限状態スペクトルモデルを提案し、音楽音響信号を楽器の単音ごとに分解できることを示した（参考文献２（M. Nakano et al.,“ Infinite-state spectrum model for music signal analysis, ” in Proc. ICASSP 2011）参照）。同様の枠組みを音響イベント検出のために導入する。
音響要素ｄの基底スペクトルの状態系列{Ｚｄ,１，．．．，Ｚｄ,Ｔ}はそれぞれ離散的な値１，．．．，Ｋｄ（すなわち状態のインデックス）をとる。これらの状態系列の同時分布は
【００３１】
【数７】
2014164126.tif 000009
【００３２】
を用いて、
【００３３】
【数８】
2014164126.tif 000010
【００３４】
となる。
【００３５】
【数９】
2014164126.tif 000011
【００３６】
であり、δ（・）はクロネッカーのデルタ関数である。ここで、π＾ｄの事前分布として、多項分布の共役事前分布である対称なディリクレ分布
【００３７】
【数１０】
2014164126.tif 000012
【００３８】
を考える。ただし、θβ（ｄ）は正のパラメータとする。このように事前分布を定めると、θβ（ｄ）に対する状態系列の同時分布は
【００３９】
【数１１】
2014164126.tif 000013
【００４０】
と書ける。このとき、｛Ｚｄ，１,...,Ｚｄ，ｔ−１，Ｚｄ，ｔ＋１,...,Ｚｄ，Ｔ｝ (以降ではＺ＾ｄ，＼ｔと表現する) が与えられたときのＺｄ，ｔの条件付き確率は次のように与えられる。
【００４１】
【数１２】
2014164126.tif 000014
【００４２】
ただし、ｎｄ，＼ｔ(ｋ)はＺｄ，ｔ’＝ｋ（Ｚｄ，ｔ’∈Ｚ＾ｄ，＼ｔ）を満たすｔ′の個数を表す。さらに、音響要素を表現するスペクトルの状態数をデータから学習する枠組みを導入する。ＣＲＰ（参考文献６（D. J. Aldous,“ Representations for partially exchangeable arrays of random variables, ”Journal of Multivariate Analysis, vol. 11, pp. 581−598, 1981.）や参考文献７（Y.W.Teh et al., “ Hierarchical Bayesian nonparametric models with applications, ”Cambridge University Press, 2010.）を参照）では、有限として考えた状態数をＫｄ→∞とすることで、Ｚｄ，ｔの条件付き確率を次のように表す。
【００４３】
【数１３】
2014164126.tif 000015
【００４４】
ここで、Ｋｔ，＼＋はｎｄ，＼ｔ(ｋ)＞０となるクラスの個数である。上式から分かるように、無限状態スペクトルモデルにおいて、各時刻に用いられる状態に着目すると、他の時刻に多く用いられている状態ほど使われやすくなる性質がある．また、新しい状態が用いられやすくなるか否かはパラメータθβ（ｄ）に影響する。
【００４５】
[１．３  ベイズロジスティック回帰に基づくマルチラベリング]
アクティベーション行列から、各フレームで発音する音響イベントを推定する（マルチラベリングする）ために、ベイズロジスティック回帰（参考文献８（D. J. C. MacKay,“ The evidence framework applied to classification networks, ”Neural Computation, vol. 4,pp. 448−472, 1992.）を参照）を利用する。これはトピックモデルを教師あり学習できるよう拡張したsLDA（参考文献９（D. M. Blei and J. D. McAuliffe,“ Supervised topic models, ”in Proc. NIPS 2007.）を参照）と類似した適用方法である。まず、フレームt における音響要素のアクティベーションをまとめてＵｔ＾＝［Ｕ１，ｔ，Ｕ２，ｔ，．．．，ＵＤ，ｔ］Ｔと表現する。また、音響イベントl に対するラベルをＸｌＮ＾＝｛Ｘｌ，ｔ１，Ｘｌ，ｔ２，．．．，Ｘｌ，ｔＮ｝（Ｘｌ，ｔｎ∈｛０，１｝） と表現すると、その尤度関数は、
【００４６】
【数１４】
2014164126.tif 000016
【００４７】
と書ける。ここで、ｗｌ＾＝［ｗｌ，１，ｗｌ，２，．．．，ｗｌ，Ｄ］Ｔは重みベクトルであり、σ（・） はロジスティックシグモイド関数
【００４８】
【数１５】
2014164126.tif 000017
【００４９】
である。さらに、次の形の単純な等方ガウス事前分布
【００５０】
【数１６】
2014164126.tif 000018
【００５１】
とαｌに対する共役超事前分布をガンマ分布
【００５２】
【数１７】
2014164126.tif 000019
【００５３】
を導入する。これは超パラメータθα（０）とθα（１）で定まる。Ｉ＾ＤはＤ×Ｄの単位行列を表す。このとき、アクティベーションベクトルＵｔ＾が与えられた下で、音響イベントラベルＸｌ，ｔの予測分布は、
【００５４】
【数１８】
2014164126.tif 000020
【００５５】
となる。上式では、ロジスティックシグモイド関数とガウス分布のたたみ込み積分を近似するために、プロビット関数に基づく近似予測分布を利用する。ｗｌＭＡＰ＾はｗ＾のＭＡＰ（最大事後確率）解を表し、実際はラベル付けされた少量の学習データを利用して推定される。
【００５６】
音響イベント書き起こしモデルのグラフィカル表現を図3 に示す。これは音響特徴行列Ｙ＾とラベル行列Ｘ＾を共に生成するモデルある。ここで、計算の都合上、基底スペクトル集合Ｈ＾の事前分布として、ガンマ分布Ｈω，ｄ（ｋ）〜Ｇａｍｍａ（θφ（ω，d），θψ（ω，d））を仮定する。また、θａ，θｂ（０），θｂ（１），θα（０），θα（１），θβ（ｄ），θφ（ω，d），θψ（ω，d）はすべて固定の超パラメータである。
また、音響特徴量行列Ｙ＾と、ある特定のフレームだけにラベルが付与された音響イベントラベル行列Ｘ＾＝[Ｘ１Ｎ＾，Ｘ２Ｎ＾，．．．，ＸＬＮ＾]Ｔとの対数尤度関数は、
【００５７】
【数１９】
2014164126.tif 000021
【００５８】
となる。
【００５９】
[２  パラメータの推論]
ｍＩＢＰのstick breaking constructionは、スライスサンプリングと動的計画法を組み合わせることで、モデルパラメータを推論できる（参考文献３（J. V. Gael et al.,“The infinite factorial hidden Markov model, ”in Proc. NIPS 2008）を参照）。スライスサンプラーが適応的に、無限次元モデルの打ち切りを行い、その後、動的計画法によって各パラメータが推論される。しかしながら、本実施の形態はロジスティック回帰に基づくマルチラベリング処理を含むため（計算コストを軽くするため）、打ち切りされたstick breaking construction の下で（Ｄを大きな値に固定して）、ギブスサンプリングによってパラメータを推論する。以下に各パラメータの推論方法を示す。
【００６０】
[２．１  Ｃω，ｔ，１，．．．，Ｃω，ｔ，Ｄの推論]
Ｃω，ｔ，１，．．．，Ｃω，ｔ，Ｄの条件付き確率は、
【００６１】
【数２０】
2014164126.tif 000022
【００６２】
となる。ここで、ｒｄ，ｔ（ｋ）＝δ（Ｚｄ，ｔ−ｋ） とし、ｐω，ｔ，ｄは
【００６３】
【数２１】
2014164126.tif 000023
【００６４】
とする。この多項分布からＣω，ｔ，１，．．．，Ｃω，ｔ，Ｄをサンプリングする。
【００６５】
[２．２  Ｚｄ，ｔの推論]
Ｚｄ，ｔの条件付き確率は、Ｙ＾の尤度関数を考慮して、
【００６６】
【数２２】
2014164126.tif 000024
【００６７】
と書ける。ここで、
【００６８】
【数２３】
2014164126.tif 000025
【００６９】
とする。これらの確率を用いて、Ｚｄ，ｔの状態割り当てを行う。
【００７０】
[２．３  Ｈω，ｄ（ｋ）の推論]
Ｈω，ｄ（ｋ）の条件付き確率は、
【００７１】
【数２４】
2014164126.tif 000026
【００７２】
となる。このガンマ分布からＨω，ｄ（ｋ）をサンプリングする。
【００７３】
[２．４  Ｕｄ，１，，．．．，Ｕｄ，Ｔの推論]
Forward-filtering backward-sampling アルゴリズム（参考文献１０（S. L. Scott, “Bayesian methods for hidden Markov models: Recursive computing in the 21st century, ”JASA, vol. 97, pp. 337−351, 2002.）を参照）を利用して、Ｕｄ，１，，．．．，Ｕｄ，Ｔを推論する。このとき、Ｕ＾のｄ行以外の値はすべて固定する。まず、ｔ＝１，．．．，Ｔに対して、
【００７４】
【数２５】
2014164126.tif 000027
【００７５】
を再帰的に計算する。ここで、Ｘｌ，ｔの尤度関数は、ラベルが付与されているフレームｔだけ計算される。次に、
【００７６】
【数２６】
2014164126.tif 000028
【００７７】
からＵｄ，Ｔをサンプリングする。そして、ｔ＝Ｔ−１，．．．，１に対して、Ｕｄ，ｔ＋１が与えられた下で、
【００７８】
【数２７】
2014164126.tif 000029
【００７９】
に従って、Ｕｄ，ｔを後方から順番にサンプリングすれば、Ｕｄ，１，，．．．，Ｕｄ，Ｔが求まる。
【００８０】
[２．５  ａ（ｄ）の推論]
ｄ＝１，．．．，Ｄ に対して、ａ（ｄ） の条件付き確率は、
【００８１】
【数２８】
2014164126.tif 000030
【００８２】
となり、このベータ分布からａ（ｄ）をサンプリングする。
【００８３】
[２．６  ｂ（ｄ）の推論]
ｄ＝１，．．．，Ｄに対して、ｂ（ｄ）の条件付き確率は、
【００８４】
【数２９】
2014164126.tif 000031
【００８５】
となり、このベータ分布からｂ（ｄ）をサンプリングする。
【００８６】
[２．７  ｗｌ＾の推論]
音響信号の一部のＮフレームには音響イベントのｏｎ／ｏｆｆに関するラベルが付いているとする。ラベルが付けられたフレームのインデックスをｔ１，ｔ２，．．．，ｔＮとし、これらのフレームの音響要素のアクティベーション｛Ｕｔ１＾，Ｕｔ２＾，．．．，ＵｔＮ＾｝ と音響イベントｌに対するラベルＸｌＮ＾＝［Ｘｌ，ｔ１，Ｘｌ，ｔ２，．．．，Ｘｌ，ｔＮ］Ｔを用いて、ｗｌを推論する。ロジスティックシグモイド関数の積である尤度関数と事前確率分布が共役な関係ではないため、ここでは、ニュートン-ラフソン法を利用してＭＡＰ（最大事後確率）推定値ｗｌＭＡＰ＾を求める。更新式は、
【００８７】
【数３０】
2014164126.tif 000032
【００８８】
と導出される。ここで、
【００８９】
【数３１】
2014164126.tif 000033
【００９０】
とする。ｗｌ（ｎｅｗ）＾とｗｌ（ｏｌｄ）＾の差が10−8 以下になったとき、ｗｌＭＡＰ＾＝ｗｌ（ｎｅｗ）＾とする。
【００９１】
[２．８  αｌの推論]
ｌ＝１，．．．，Ｌに対して、αｌの条件付き確率は、
【００９２】
【数３２】
2014164126.tif 000034
【００９３】
と近似され、このガンマ分布からαｌをサンプリングする。
【００９４】
＜音響信号分析装置のシステム構成＞
次に、音響信号の信号パラメータを分析して出力する音響信号分析装置に、本発明を適用した場合を例にして、本発明の実施の形態を説明する。
【００９５】
図４に示すように、本実施の形態に係る音響信号分析装置は、ＣＰＵと、ＲＡＭと、後述する音響信号分析処理ルーチンを実行するためのプログラムを記憶したＲＯＭとを備えたコンピュータで構成され、機能的には次に示すように構成されている。
【００９６】
本実施の形態に係る音響信号分析装置は、入力部１と、記憶部２と、演算部３と、出力部４とを備えている。本実施の形態は、分析対象である音響信号の時系列の一部分に対してＬ個の音響イベントの値をラベル付けしたデータを学習データとして入力する。そして、前述したパラメータ推定アルゴリズムを用いて信号解析を行い、分析対象である音響信号の時系列の全ての部分に対して、推定したパラメータを適用する構成である。
【００９７】
入力部１には、分析対象である音響信号の時系列と、当該音響信号の時系列について、ある特定のフレームだけにラベルが付与された音響イベントラベル行列Ｘ＾＝[Ｘ１Ｎ＾，Ｘ２Ｎ＾，．．．，ＸＬＮ＾]Ｔとが入力される。ここで、Ｘ＾は、音響信号の時系列の一部分（ある特定のフレーム）について予め付与された教師ラベル行列である。
【００９８】
記憶部２は、入力部１に入力された音響信号の時系列と、Ｘ＾とを記憶する。また、記憶部２には、予め定められた超パラメータθａ＝１，θｂ（０）＝１０，θｂ（１）＝１，θβ（ｄ）＝１，θφ（ω，d）＝１，θψ（ω，d）＝１  （ｄ ＝１，．．．，Ｄ，ω＝１，．．．，Ω），θα（０）＝１，θα（１）＝１０００が記憶されている。これらの超パラメータの各値は予備実験によって決定した。また、記憶部２は、後述する各処理での結果等を記憶する。
【００９９】
演算部３は、音響特徴抽出部３０と、パラメータ初期値生成部３２と、モデルパラメータ更新部３４と、音響イベントラベル推論部３６とから構成されている。
【０１００】
音響特徴抽出部３０は、記憶部２に記憶された音響信号の時系列を読み出し、当該音響信号に基づいて、音響特徴量行列Ｙ＾（以下、Ｙ＾と称する。Ｙ＾は二次元配列である。）を作成する。本実施の形態では、音響イベントの音色に着目し、フレーム長１００ｍｓ、フレームシフト長１００ｍｓとして、フレームごとに短時間フーリエ変換を行い、その振幅スペクトルをメルフィルタバンク処理して得られる出力値であるＹω,tを各要素にもつ二次元配列Ｙ＾（ω＝１,...,Ω、ｔ＝１,...,Ｔ）を生成し、記憶部２へ格納する。このとき、例えば、メルフィルタの数を２４個とすれば、Ω＝２４であり、Ｔはフレームの総数に相当する。また、上記（１）式に示すように、本実施の形態では、ポアソン分布に従う確率変数として、音響特徴量行列をモデル化するため、Ｙ＾の全ての要素の値を整数値に丸め込む。
【０１０１】
パラメータ初期値生成部３２は、基底スペクトル集合Ｈ＾（以下、Ｈ＾と称する。）、アクティベーション集合Ｕ＾（以下、Ｕ＾と称する。）、状態系列Ｚ＾（以下、Ｚ＾と称する。）、重みベクトルｗｌ＾（以下、ｗｌ＾と称する。）、精度パラメータαｌ、及び遷移確率ａ＾並びにｂ＾の各パラメータの初期値を生成し、記憶部２へ格納する。
なお、Ｈ＾は三次元配列、Ｕ＾とＺ＾とは二次元配列、ｗｌ＾とａ＾とｂ＾とは一次元配列である。
【０１０２】
パラメータ初期値生成部３２は、音響特徴抽出部３０により作成されたＹ＾に基づいて、上記（１）式に従って、通常のＮＭＦを適用し、推定されたＨ＾とＵ＾を初期値とし、記憶部２に格納する。ここで、推定されたＨ１：Ω，ｄはｄ番目の音響要素の初期基底スペクトルとする（初期状態数はＫｄ＝１とする）。
【０１０３】
具体的には、パラメータ初期値生成部３２は、Ｈ＾について、各音響要素ｄの基底の状態ｋを１のみとする。そして、Ｈ＾の要素である全てのパラメータＨω,d(1)が非負値であり、Ｕ＾の要素である全てのパラメータＵd,tが非負値である、という条件の下で、パラメータＨω,d(1)とパラメータＵd,tとの積を、全てのｄについて足し合わせたモデルについて、Ｙ＾とモデルとの距離を表わした目的関数の値が小さくなるように、パラメータＨω,d(1)及びパラメータＵd,tの各々を推定して、Ｈ＾及びＵ＾の初期値を設定する。
【０１０４】
ＮＭＦは、周知技術により実現できる（例えば、参考文献１１（A.T.Cemgil,“Bayesian inference in non-negative matrix factorisation models, ”in University of Cambridge, 2008.） や、参考文献１２（M. Hoffman, D. Blei, and P. Cook,“ Bayesian nonparametric matrix factorization for recorded music, ”in Proc. ICML,2010.）を参照）ため、詳細な説明を省略する。
【０１０５】
また、Ｕ＾については、その中央値よりも大きい要素は１に、中央値よりも小さい要素は０に二値化して初期値とする。
【０１０６】
Ｚ＾については、全ての要素の初期値として１が設定され、記憶部２に格納される。
【０１０７】
音響要素数Ｄについては、適切な初期値が設定され、記憶部２に格納される。本実施の形態では、打ち切りstick breaking construcitonをＵ＾に適用するため、音響要素数Ｄの値はできる限り大きい値が望ましい。後述する評価実験では、Ｄ＝１００とした。
【０１０８】
遷移確率ａ＾については、上記（６）式に従って、ａ＾＝{ａ１,...,ａＤ}についての初期値が生成され、そして遷移確率ｂ＾については、上記（５）式に従って、初期値が生成され、各々の初期値が記憶部２に格納される。
【０１０９】
ｗｌ＾、精度パラメータαｌについては、音響イベントｌ＝１，．．．，Ｌに対して、初期値として適切な値が設定され、各々の初期値が記憶部２に格納される。
【０１１０】
モデルパラメータ更新部３４は、前述したパラメータ推定アルゴリズムを用いて、各パラメータを更新して各パラメータを推定する。モデルパラメータ更新部３４は、音源分離推論部３４０と、状態系列推論部３４１と、基底スペクトル推論部３４２と、アクティベーション推論部３４３と、遷移確率推論部３４４と、重みベクトル推論部３４５と、精度パラメータ推論部３４６と、収束判定部３４７とから構成されている。
【０１１１】
音源分離推論部３４０は、上記（１８）式及び（１９）式に従って、数値的サンプリングによって、各音響要素ｄについての音響特徴量成分を表す行列Ｃ＾（以下、Ｃ＾と称する。Ｃ＾は三次元配列である。）のすべての要素を推論する。具体的には、記憶部２に記憶されたＹ＾、Ｚ＾、Ｈ＾、及びＵ＾に基づいて、上記（１８）式及び（１９）式に従って、Ｃ＾のすべての要素をサンプリングして推論し、記憶部２に格納する。
【０１１２】
状態系列推論部３４１は、上記（２０）式及び（ｃ）式から、ギブスサンプリングによってＺ＾のすべての要素の推論を行う。具体的には、記憶部２に記憶されたＹ＾、超パラメータθφ（ω，d），θψ（ω，d）（ｄ ＝１，．．．，Ｄ，ω＝１，．．．，Ω）、Ｃ＾、Ｚ＾、及びＵ＾に基づいて、上記（２０）式及び（ｃ）式に従って、ｐ（Ｚｄ,１,...,Ｚｄ,Ｔ｜Ｃ・,１,ｄ,...,Ｃ・,Ｔ,ｄ）の値の大きなＺｄ,１,...,Ｚｄ,Ｔを確率的に求め、記憶部２に格納する。なお、Ｚｄ,ｔを求める際のＺ＾については、Ｚｄ,ｔ以外の｛Ｚｄ,１,...,Ｚｄ,ｔ−１,Ｚｄ,ｔ＋１,...,Ｚｄ,ｔ｝（＝Ｚｄ,＼ｔ）の値を既知として、Ｚｄ,ｔをサンプリングして推論し、記憶部２に格納する。
【０１１３】
基底スペクトル推論部３４２は、上記（２１）式から、数値的サンプリングによってＨ＾のすべての要素を推論する。具体的には、記憶部２に記憶された超パラメータθφ（ω，d），θψ（ω，d）（ｄ ＝１，．．．，Ｄ，ω＝１，．．．，Ω）、Ｃ＾、Ｚ＾、及びＵ＾に基づいて、上記（２１）式に従って、数値的サンプリングによってＨ＾のすべての要素をサンプリングして推論し、記憶部２に格納する。
【０１１４】
アクティベーション推論部３４３は、上記（２２）式及び（２３）式からなるForward filtering backward sampling アルゴリズムに基づいて、Ｕ＾のすべての要素を推論する。具体的には、記憶部２に記憶されたＹ＾、Ｃ＾、Ｚ＾、Ｈ＾、Ｘ＾及びＵ＾に基づいて、まずForward filteringによって、上記（２２）式に従って、各ｔに対するＵｄ,ｔの事後分布を再帰的に計算した後、Ｂackward samplingによって、上記（２３）式に従って、各ｔに対するＵｄ,ｔを数値的にサンプリングして求め、記憶部２に格納する。なお、Forward filteringでは、Ｕｄ,ｔの事後分布を求める際に、Ｕｄ,ｔ-1の事後分布を用い、Ｂackward samplingでは、Ｕｄ,ｔを求める際に、Ｕｄ,ｔ＋1を用いる。
【０１１５】
遷移確率推論部３４４は、記憶部２に記憶された超パラメータθｂ（０）、θｂ（１）及びＵ＾に基づいて、上記（２４）式及び（２５）式に従って、遷移確率ａ（ｄ）及びｂ（ｄ）を推論する。具体的には、記憶部２に記憶された遷移確率ａ＾、及びＵ＾に基づいて、上記（２４）式に従って、ｄ＝１,...,Ｄに対するａ（ｄ）をサンプリングして推論し、記憶部２に格納する。また、ｂ（ｄ）については、記憶部２に記憶された超パラメータθｂ（０）、θｂ（１）、及びＵ＾に基づいて、上記（２５）式に従って、ｄ＝１,...,Ｄに対するｂ（ｄ）をサンプリングにより推論し、記憶部２に格納する。
【０１１６】
重みベクトル推論部３４５は、音響イベントｌ＝１，．．．，Ｌに対して、ニュートン-ラフソン法を利用して、最大事後確率推定値ｗｌＭＡＰ＾を求める。具体的には、音響イベントｌ＝１，．．．，Ｌに対して、記憶部２に記憶されたＵ＾、ＸｌＮ＾、精度パラメータαｌ、及びｗｌ（ｏｌｄ）＾に基づいて、上記（２６）式に従って、重みベクトルｗｌ（ｎｅｗ）＾を更新することを繰り返し、更新したｗｌ（ｎｅｗ）＾とｗｌ（ｏｌｄ）＾の差が１０ー８以下になったとき、ｗｌＭＡＰ＾＝ｗｌ（ｎｅｗ）＾として最大事後確率推定値ｗｌＭＡＰ＾を更新し、記憶部２に格納する。
【０１１７】
精度パラメータ推論部３４６は、音響イベントｌ＝１，．．．，Ｌに対して、記憶部２に記憶された超パラメータθα（０）、θα、及びｗｌＭＡＰ＾に基づいて、上記（２８）に従って、数値的サンプリングによって、αｌを推論し、記憶部２に格納する。
【０１１８】
収束判定部３４７は、記憶部２に記憶されたＹ＾及びＸ＾と、推論されたＺ＾、Ｈ＾、Ｕ＾、ｗ＾とを用いて、上記（１７）式の対数尤度関数を計算する。更新前のモデルパラメータを用いて計算した対数尤度関数の値と更新後の値との誤差が、所定の閾値以下であれば、モデルパラメータの更新が収束したと判定する。一方、誤差が所定の閾値を超えるようであれば、モデルパラメータの更新は収束していないと判定し、音源分離推論部３４０の処理に戻る。本実施の形態ではこの誤差の閾値を実験的にε＝１．０×１０−５とした。
【０１１９】
音響イベントラベル推論部３６は、Ｌ個の音響イベントｌの各々について、記憶部２に記憶されたＵ＾、ＸｌＮ^、精度パラメータαｌ及びｗｌＭＡＰ＾に基づいて、上記（１６）式に従って、各時刻ｔのパラメータＸｌ,tの予測確率を算出し、時刻ｔ毎に、時刻ｔの予測確率が、０．５より大きい場合には、パラメータＸｌ,tの値を１とし、時刻ｔの予測確率が、０．５以下である場合には、パラメータＸｌ,tの値を０とする。
【０１２０】
出力部４は、記憶部２に記憶されているモデルパラメータＣ＾、Ｚ＾、Ｈ＾、Ｕ＾、ａ＾、ｂ＾、ｗｌＭＡＰ＾をすべて出力すると共に、音響イベントラベル推論部３６により求められた各音響イベントｌに対する時刻ｔ毎のパラメータＸｌ,tを出力する。
【０１２１】
なお、収束したか否かを判定する方法としては、対数尤度関数を用いる方法以外に、モデルパラメータ各々の値を更新前と更新後とで比較しても良いし、予め定めた繰り返し回数に到達したか否かで判定を行っても良い。本実施の形態ではモデルパラメータ各々の値を更新前と更新後とで比較する場合、この誤差の閾値がε＝１．０×１０−５であれば良好な結果であることを実験的に確認した。また、予め定めた繰り返し回数に到達したか否かで判定をする場合、１０００回の繰り返し回数が必要であることも実験的に確認している。
【０１２２】
＜音響信号分析装置の作用＞
次に、本実施の形態に係る音響信号分析装置の作用について説明する。まず、分析対象の時系列信号として音響信号が音響信号分析装置に入力され、記憶部２に格納される。また、音響イベントラベル行列Ｘ＾が音響信号分析装置に入力され、記憶部２に格納される。そして、音響信号分析装置において、図５に示す音響信号分析処理ルーチンが実行される。
【０１２３】
まず、ステップＳ１００において、記憶部２から、各フレーム内の音響信号の時系列を読み込む。そして、音響信号の時系列に対して、短時間フーリエ変換を用いた時間周波数分析を行い、振幅スペクトルを結果として得る。そして、その振幅スペクトルに対し、メルフィルタバンク処理して得られる音響特徴量Ｙω,tを各（ω，ｔ）の要素にもつ二次元配列Ｙ＾を生成し、記憶部２へ格納する。
【０１２４】
そして、ステップＳ１０２において、記憶部２に記憶されている超パラメータθａ、θｂ（０）、θｂ（１）、θβ（ｄ）、θφ（ω，d）、θψ（ω，d）（ｄ ＝１，．．．，Ｄ，ω＝１，．．．，Ω）、θα（０）、θα（１）の値を読み込む。
【０１２５】
次にステップＳ１０６において、パラメータ初期値生成部３２によって、Ｈ＾、Ｕ＾、Ｚ＾、音響要素数Ｄ、遷移確率ａ＾、ｂ＾、ｗｌ＾、及び精度パラメータαｌの各パラメータの初期値を生成する。
具体的には、上記（１）式に従って、上記ステップＳ１００において生成されたＹ＾に基づいて、Ｈ＾とＵ＾を推定し、初期値として記憶部２へ記憶する。Ｚ＾については、全ての要素の初期値として１を設定し、音響要素数Ｄについては、初期値として例えばＤ＝１００を設定する。また遷移確率ａ＾については、上記（６）式に従って、ａ＾＝{ａ１,...,ａＤ}についての初期値が生成され、遷移確率ｂ＾については、上記（５）式に従って初期値が生成される。ｗｌ＾、精度パラメータαｌについては、音響イベントｌ＝１，．．．，Ｌに対して、初期値として適切な値が設定される。そして、生成された各々の初期値を記憶部２へ記憶する。
【０１２６】
ステップＳ１０８では、音源分離推論部３４０によって、Ｃ＾の全ての要素を推論する。具体的には、上記ステップＳ１００において生成されたＹ＾、並びに上記ステップＳ１０６で設定されたＺ＾、Ｈ＾及びＵ＾、又は、後述するステップＳ１１０〜ステップＳ１１４で前回推論されたＺ＾、Ｈ＾及びＵ＾に基づいて、上記（１８）式及び（１９）式に従って、Ｃ＾のすべての要素をサンプリングして推論し、記憶部２へ記憶する。
【０１２７】
ステップＳ１１０では、状態系列推論部３４１によって、ギブスサンプリングによって、Ｚ＾のすべての要素の推論を行う。具体的には、上記ステップＳ１００で生成されたＹ＾、上記ステップＳ１０２で読み込まれた超パラメータθφ（ω，d），θψ（ω，d）（ｄ ＝１，．．．，Ｄ，ω＝１，．．．，Ω）、上記ステップＳ１０６で設定されたＺ＾及びＵ＾、又は本ステップＳ１１０及び後述するステップＳ１１４で前回推論されたＺ＾及びＵ＾、並びに上記ステップＳ１０８で推論されたＣ＾に基づいて、上記（２０）式及び（ｃ）式に従って、Ｚｄ,１,...,Ｚｄ,Ｔをサンプリングして推論する。
なお、Ｚｄ,ｔを求める際の状態系列Ｚ＾については、Ｚｄ,ｔ以外の｛Ｚｄ,１,...,Ｚｄ,ｔ−１,Ｚｄ,ｔ＋１,...,Ｚｄ,ｔ｝（＝Ｚｄ,＼ｔ）の値を既知として、Ｚｄ,ｔをサンプリングして推論する。ここで、Ｚｄ,１,...,Ｚｄ,ｔ−１までは、本ステップＳ１１０で推論された直前の値を用い、Ｚｄ,ｔ＋１,...,Ｚｄ,ｔについては、上記ステップＳ１０６で生成されたＺ＾又は本ステップＳ１１０で前回推論されたＺ＾を用いて、Ｚｄ,ｔをサンプリングして推論し、記憶部２へ記憶する。
【０１２８】
ステップＳ１１２では、基底スペクトル推論部３４２によって、Ｈ＾の全ての要素の推論を行う。具体的には、上記ステップＳ１０２で読み込まれた超パラメータθφ（ω，d），θψ（ω，d）（ｄ ＝１，．．．，Ｄ，ω＝１，．．．，Ω）、上記ステップＳ１０６で設定されたＵ＾、又は後述するステップＳ１１４で前回推論されたＵ＾、上記ステップＳ１０８で推論されたＣ＾、及び上記ステップＳ１１０で推論されたＺ＾に基づいて、上記（２１）式に従って、Ｈ＾のすべての要素をサンプリングして推論し、記憶部２へ記憶する。
【０１２９】
ステップＳ１１４では、アクティベーション推論部３４３によって、Ｕ＾の全ての要素の推論を行う。具体的には、上記ステップＳ１００において生成されたＹ＾、記憶部２に記憶されたＸｌＮ＾、上記ステップＳ１０８で推論されたＣ＾、上記ステップＳ１１０で推論されたＺ＾、上記ステップＳ１１２で推論されたＨ＾、及び上記ステップＳ１０６で設定されたＵ＾又は本ステップＳ１１４で前回推論されたＵ＾に基づいて、上記（２２）式及び（２３）式に従って、Ｕ＾の全ての要素の推論を行う。詳細には、まずForward filteringによって、上記（２２）式に従って、各ｔに対するＵｄ,ｔの事後分布を再帰的に計算した後、Ｂackward samplingによって、上記（２３）式に従ってＵｄ,ｔを数値的にサンプリングして求めることにより、Ｕ＾の全ての要素を推論し、記憶部２へ記憶する。
【０１３０】
ステップＳ１１６では、遷移確率推論部３４４によって、遷移確率ａ（ｄ）、ｂ（ｄ）の推論を行う。具体的には、本ステップＳ１１６で今回または前回生成された遷移確率ａ、上記ステップＳ１０２で読み込まれた超パラメータθｂ（０）、θｂ（１）、及び上記ステップＳ１１４で推論されたＵｄ,１,..., Ｕｄ,Ｔに基づいて、上記（２４）式、及び（２５）式に従って、遷移確率ａ（ｄ）、ｂ（ｄ）の推論を行う。まず、上記（２４）式に基づいて、ｄ＝１,...,Ｄに対するａ（ｄ）を推論する。ｂ（ｄ）については、上記（２５）式に基づいて、ｄ＝１,...,Ｄに対するｂ（ｄ）を推論する。そして、推論された遷移確率ａ（ｄ）、ｂ（ｄ）を記憶部２へ記憶する。
【０１３１】
ステップＳ１１８において、重みベクトル推論部３４５によって、ｗｌＭＡＰ＾の推論を行う。具体的には、音響イベントｌ＝１，．．．，Ｌに対して、上記ステップＳ１１４で推論されたＵ＾、記憶部２に記憶されたＸｌＮ＾、上記ステップＳ１０６で設定されたか、又は後述するステップＳ１２０で前回推論された精度パラメータαｌ、及び上記ステップＳ１０６で設定されたか、又は本ステップＳ１１８で前回推論されたｗｌ（ｏｌｄ）＾に基づいて、上記（２６）式に従って、重みベクトルｗｌ（ｎｅｗ）＾を更新することを繰り返し、更新したｗｌ（ｎｅｗ）＾とｗｌ（ｏｌｄ）＾の差が１０ー８以下になったとき、ｗｌＭＡＰ＾＝ｗｌ（ｎｅｗ）＾として最大事後確率推定値ｗｌＭＡＰ＾を更新し、記憶部２に格納する。
【０１３２】
ステップＳ１２０において、精度パラメータ推論部３４６によって、αｌを推論する。
具体的には、音響イベントｌ＝１，．．．，Ｌに対して、上記ステップＳ１０２で読み込まれた超パラメータθα（０）、θα（１）、及び上記ステップＳ１１８で推論されたｗｌＭＡＰ＾に基づいて、上記（２８）に従って、数値的サンプリングによって、αｌを推論する。そして、推論されたαｌを記憶部２へ記憶する。
【０１３３】
ステップＳ１２２では、収束判定部３４７によって、上記ステップＳ１００で得られたＹ＾及び記憶部２に記憶されたＸ＾と、上記ステップＳ１１０〜Ｓ１１８で推論されたＺ＾、Ｈ＾、Ｕ＾、及びｗｌ＾とを用いて、上記（１７）式の対数尤度関数を計算する。そして、更新前のモデルパラメータを用いて計算した対数尤度関数の値と更新後の値との誤差が、所定の閾値以下であれば、収束したと判定する。収束していないと判定された場合には、ステップＳ１０８へ移行し、上記ステップＳ１０８〜ステップＳ１２０で推論したパラメータを用いて、上記ステップＳ１０８〜ステップＳ１２０の処理を繰り返す。収束したと判定された場合には、ステップＳ１２４へ移行する。
【０１３４】
ステップＳ１２４において、音響イベントラベル推論部３６によって、Ｌ個の音響イベントｌの各々について、記憶部２に記憶されたＵ＾、ＸｌＮ^、精度パラメータαｌ及びｗｌＭＡＰ＾に基づいて、上記（１６）式に従って、各時刻ｔのパラメータＸｌ,tの予測確率を算出する。そして、Ｌ個の音響イベントｌの各々について、時刻ｔ毎に、時刻ｔの予測確率が、０．５より大きい場合には、パラメータＸｌ,tの値を１とし、時刻ｔの予測確率が、０．５以下である場合には、パラメータＸｌ,tの値を０とする。
【０１３５】
ステップＳ１２６では、出力部４によって、推論された各パラメータ（Ｃ＾、Ｚ＾、Ｈ＾、Ｕ＾、ａ＾、ｂ＾、ｗｌＭＡＰ＾）を結果として出力すると共に、音響イベントラベル推論部３６により求められた各音響イベントｌの時刻ｔ毎のパラメータＸｌ,tを出力し、音響信号分析処理ルーチンが終了する。
【０１３６】
＜実験結果＞
英語学習用ポッドキャストの音響信号（計３５分のうち、はじめの５分）を用いて、本技術の性能を評価する。音響信号は１６ｋＨｚサンプリングで１６ビット量子化されたものである。
【０１３７】
図６の上段の図は、１００ｍｓごとに手動でラベル付けされた音響イベントラベルを示す。黒色が音響イベントのｏｎを、白色がｏｆｆを表す。この音響信号は“音楽”、“効果音”、“電話のベル音”、５名の“音声”からなる８種類の音響イベントを含む。音響信号は、フレームシフト長１００ｍｓ、フレーム長１００ｍｓ、ハニング窓を用いてフレームに分割され、短時間フーリエ変換によって振幅スペクトログラムに変換される。そして、各フレームの振幅スペクトルをフィルタバンク処理して得られる出力値（２４個の値）を、音響特徴量行列とする（すなわち、Ｙ＾は２４×３０００の行列である）。
【０１３８】
図６の上段の図に示されるように、学習データ（あらかじめラベル付けされたデータ）として３ 種類（（１）５０秒、（２）１００秒、（３）１５０秒）を用意した。そして、後半の１５０秒に対する推定結果に基づいて、本技術の性能を評価する。評価尺度として、各イベントのラベル予測確率｛ｐ（Ｘ｛ｌ，１５０１｝｜Ｕ＾ｔ，Ｘ＾ｌＮ）， ｐ（Ｘ｛ｌ，１５０２｝｜Ｕ＾ｔ，Ｘ＾ｌＮ），…ｐ（Ｘ｛ｌ，３０００｝｜Ｕ＾ｔ,Ｘ＾ｌＮ）｝と正解ラベルとを比較することで計算されるＲＯＣ曲線の下側面積ＡＵＣの値を利用する。
【０１３９】
図６の中段及び下段の図は、アクティベーション行列Ｕ＾と学習データ（３）を利用して推定されたラベル付け結果を示す。
Ｘｌ，ｔ＝Ｉ（ｐ（Ｘ｛ｌ,ｔ｝|Ｕ＾ｔ，Ｘ＾ｌＮ)＞０．５）
の基準に基づいて、各フレームに付与する音響イベントｌのラベルの値が０であるか１であるかを判定した。ここで、ｐ（Ｘ｛ｌ,ｔ｝｜Ｕ＾ｔ，Ｘ＾ｌＮ)は予測確率である。また、Ｉ（Ａ）はインジケータ関数であり、条件Ａが真であれば１を、偽であれば０を返す関数である。アクティベーション行列のロジスティック回帰によって推定されるラベル付結果が、正解ラベルと全体的に近い結果が得られたことがわかる。また、女性の話者を識別することが難しいこともわかる。
【０１４０】
図７は、学習データの長さに対するラベル付性能を評価した図である。ここで、“電話のベル音” と“男性Ａの音声” は評価データに含まれていないため、これらについては評価しない。学習データの長さが増えるにつれて、各イベントのラベル付性能を表すＡＵＣの値が向上することがわかる。
【０１４１】
図８は、ベースライン手法と本実施の形態の性能を比較する図である。どちらも学習データ（３）を使用した結果である。ベースライン手法では、音響イベントのラベルが付与されるフレームの音響特徴量（２４個のメルフィルタバンク出力値）の分布をガウス混合モデル（ＧＭＭ）を利用して学習し、評価データに対して、各音響イベントラベルの事後確率を計算する。そして、５点移動平均フィルタを利用して、事後確率を平滑化し、平均ＡＵＣの値を計算した。本技術はベースライン手法と同等以上の性能が得られていることがわかった。ただし注目すべきは、（１）や（２）のように学習データが少ない場合、ベースライン手法はＧＭＭを効果的に学習できない。しかしながら、本技術は、たとえ少量の学習データであったとしても、ラベル付けされていないデータをうまく活用することによって、全体のイベントラベルを推定できることが分かった。
【０１４２】
また、識別モデルにとって、学習データに含まれていない未知の音響イベントを検出することは難しい。しかし、生成的なアプローチである本技術が未知の音響イベントを検出できる可能性があることを示す。ここでは、“効果音”に着目する。実際、この効果音は様々な種類の音に対してラベル付けされているため、性能は低い（上記図７参照）。しかし、推定されたアクティベーション行列を見てみると興味深い（図９参照）。“効果音” のラベルが付与される区間において、特定の音響要素がアクティブとなっていることがわかる。これらの音響要素が新しい音響イベントとして検出されていることを意味する。
【０１４３】
以上説明したように、本発明の実施の形態に係る音響信号分析装置によれば、音響要素ｄの各々に対する音響特徴量成分Ｃω,t,ｄを要素にもつＣ＾をサンプリングにより決定し、音響要素ｄの各々に対して各時刻ｔに基底の状態が何れであるか、または新規の状態であることを表わすパラメータＺｄ,tを要素にもつＺ^をサンプリングにより決定し、音響要素ｄの各々における基底の各状態ｋの基底スペクトルを表わすパラメータＨω,d(k)を要素にもつＨ^をサンプリングにより決定し、予め与えられたＬ個の音響イベントｌの各々に対する一部の時刻ｔnの各々の音響イベントの有無を表すパラメータＸｌ,tｎ（＝｛０、１｝）を要素にもつＸｌＮ^を用いて、音響要素ｄの各々に対する各時刻ｔの発音の有無を表すパラメータＵd,tを要素にもつＵ^をサンプリングにより決定し、音響要素ｄの各々について、パラメータＵd,tが０から１へ遷移する遷移確率ａ（d）をサンプリングにより決定すると共にパラメータＵd,tが１から１へ遷移する遷移確率ｂ（d）をサンプリングにより決定し、予め与えられたＸｌＮ^を用いて、音響イベントｌの各々に対するＤ個の音響要素ｄの各々の重みｗｌ，ｄを要素にもつｗｌ＾を求め、音響イベントｌの各々に対する精度パラメータαｌを算出することを、予め定められた収束条件を満たすまで繰り返し行い、Ｃ＾、Ｚ^、Ｈ^、Ｕ^及びｗｌ＾を出力することにより、音響信号の時系列データから、そこに含まれる音響イベントが、どのような音響要素で構成されているのかを精度よく推定することができる。
【０１４４】
また、本発明の実施の形態に係る音響信号分析装置によれば、音響信号の時系列データ及び一部のフレームに対する音響イベントのラベル行列から、全ての音響イベントの音響的特徴とその発音区間を自律的に推定することができる。
【０１４５】
なお、本発明は、上述した実施形態に限定されるものではなく、この発明の要旨を逸脱しない範囲内で様々な変形や応用が可能である。
【０１４６】
例えば、本実施の形態における音響特徴抽出部３０では、音響信号の時系列に対し、フレームごとに短時間フーリエ変換を行い、その振幅スペクトルをメルフィルタバンク処理してＹ＾を生成する場合を例に説明したが、これに限定されるものではなく、振幅スペクトログラムそのものを音響特徴量行列Ｙ＾としてもよい。
【０１４７】
また、上述の音響信号分析装置は、内部にコンピュータシステムを有しているが、「コンピュータシステム」は、ＷＷＷシステムを利用している場合であれば、ホームページ提供環境（あるいは表示環境）も含むものとする。
【０１４８】
また、本願明細書中において、プログラムが予めインストールされている実施形態として説明したが、当該プログラムを、コンピュータ読み取り可能な記録媒体に格納して提供することも可能である。
【符号の説明】
【０１４９】
１     入力部
２     記憶部
３     演算部
４     出力部
３０   音響特徴抽出部
３２   パラメータ初期値生成部
３４   モデルパラメータ更新部
３６   音響イベントラベル推論部
３４０ 音源分離推論部
３４１ 状態系列推論部
３４２ 基底スペクトル推論部
３４３ アクティベーション推論部
３４４ 遷移確率推論部
３４５ 重みベクトル推論部
３４６ 精度パラメータ推論部
３４７ 収束判定部
【図１】
2014164126.tif 000035
【図２】
2014164126.tif 000036
【図３】
2014164126.tif 000037
【図４】
2014164126.tif 000038
【図５】
2014164126.tif 000039
【図６】
2014164126.tif 000040
【図７】
2014164126.tif 000041
【図８】
2014164126.tif 000042
【図９】
2014164126.tif 000043