(57)【特許請求の範囲】
【請求項１】
学習サンプルのソース音声の時系列データとターゲット音声の時系列データとからなるパラレルデータを入力として、前記ソース音声の時系列データから抽出される各時刻のスペクトル特徴量ベクトルと、前記ターゲット音声の時系列データから抽出される、各時刻の基本周波数とに基づいて、前記ソース音声の各時刻のスペクトル特徴量ベクトルと、前記ターゲット音声の各時刻の基本周波数との間の関係をモデル化した第１確率分布のパラメータを学習する第１モデルパラメータ学習部と、
前記ターゲット音声の各時刻の基本周波数に基づいて、基本周波数パターン生成過程をモデル化した第２確率分布のパラメータを学習する第２モデルパラメータ学習部と、
予測対象のソース音声の時系列データを入力として、前記予測対象のソース音声の時系列データから抽出される各時刻のスペクトル特徴量ベクトルと、前記第１モデルパラメータ学習部によって学習された前記第１確率分布のパラメータと、前記第２モデルパラメータ学習部によって学習された前記第２確率分布のパラメータとに基づいて、前記第１確率分布と前記第２確率分布とを用いて表される規準を大きくするように、前記予測対象のソース音声に対応する前記ターゲット音声の各時刻の基本周波数を予測する基本周波数予測部と、
を含む基本周波数パターン予測装置。
【請求項２】
前記第１確率分布を、前記ソース音声の各時刻のスペクトル特徴量ベクトルと前記ターゲット音声の各時刻の基本周波数と前記基本周波数の動的成分との同時確率分布を表す混合正規分布とし，
前記第２確率分布を、
各時刻の基本周波数と、
隠れマルコフモデルの各時刻の状態からなる状態系列、又は各時刻における甲状軟骨の平行移動運動によって生じる基本周波数パターンを表すフレーズ指令及び甲状軟骨の回転運動によって生じる基本周波数パターンを表すアクセント指令のペアからなる指令関数との組み合わせの確率分布とした請求項１記載の基本周波数パターン予測装置。
【請求項３】
前記規準を、
前記第１確率分布と前記第２確率分布との積を用いて表される、各時刻の基本周波数と、前記状態系列との組み合わせに応じた関数、または
前記第１確率分布と前記第２確率分布との積を用いて表される、各時刻の基本周波数と、各時刻の前記指令関数との組み合わせに応じた関数とした請求項２記載の基本周波数パターン予測装置。
【請求項４】
前記第１モデルパラメータ学習部は、ＥＭ（Expectation-Maximization）アルゴリズムにより、前記第１確率分布から求められる前記ソース音声の各時刻のスペクトル特徴量ベクトル及び前記ターゲット音声の各時刻の基本周波数の尤もらしさが大きくなるように、前記第１確率分布のパラメータを学習する請求項１〜請求項３の何れか１項記載の基本周波数パターン予測装置。
【請求項５】
前記第２モデルパラメータ学習部は、ＥＭ（Expectation-Maximization）アルゴリズムにより、前記第２確率分布から求められる、各時刻の基本周波数と、隠れマルコフモデルの各時刻の状態からなる状態系列との尤もらしさが大きくなるように、前記第２確率分布のパラメータとして、前記状態系列における状態遷移確率、及び各時刻における状態に応じたフレーズ指令の振幅及び各アクセント指令の振幅を表すパラメータ群を学習するか、又は
各時刻の基本周波数が与えられたときの、各時刻のフレーズ指令及びアクセント指令のペアからなる指令関数及び前記パラメータ群の対数事後確率を目的関数として、前記目的関数を増加させるように、前記指令関数及び前記パラメータ群を、前記第２確率分布のパラメータとして学習する請求項１〜請求項４の何れか１項記載の基本周波数パターン予測装置。
【請求項６】
前記基本周波数予測部は、
前記第１確率分布と前記第２確率分布との積を用いて表される、各時刻の基本周波数と、隠れマルコフモデルの各時刻の状態からなる状態系列との組み合わせに応じた関数が大きくなるように、各時刻の基本周波数と前記状態系列とを推定することにより、前記予測対象のソース音声に対応する前記ターゲット音声の各時刻の基本周波数を予測するか、又は
前記第１確率分布と前記第２確率分布との積を用いて表される、各時刻の基本周波数、各時刻のフレーズ指令及びアクセント指令のペアからなる指令関数との組み合わせに応じた関数が大きくなるように、各時刻の基本周波数と、各時刻の前記指令関数を推定することにより、前記予測対象のソース音声に対応する前記ターゲット音声の各時刻の基本周波数を予測する請求項１〜請求項５の何れか１項記載の基本周波数パターン予測装置。
【請求項７】
第１モデルパラメータ学習部と、第２モデルパラメータ学習部と、基本周波数予測部とを含む基本周波数パターン予測装置における基本周波数パターン予測方法であって、
前記第１モデルパラメータ学習部が、学習サンプルのソース音声の時系列データとターゲット音声の時系列データとからなるパラレルデータを入力として、前記ソース音声の時系列データから抽出される各時刻のスペクトル特徴量ベクトルと、前記ターゲット音声の時系列データから抽出される、各時刻の基本周波数とに基づいて、前記ソース音声の各時刻のスペクトル特徴量ベクトルと、前記ターゲット音声の各時刻の基本周波数との間の関係をモデル化した第１確率分布のパラメータを学習し、
前記第２モデルパラメータ学習部が、前記ターゲット音声の各時刻の基本周波数に基づいて、基本周波数パターン生成過程をモデル化した第２確率分布のパラメータを学習し、
前記基本周波数予測部が、予測対象のソース音声の時系列データを入力として、前記予測対象のソース音声の時系列データから抽出される各時刻のスペクトル特徴量ベクトルと、前記第１モデルパラメータ学習部によって学習された前記第１確率分布のパラメータと、前記第２モデルパラメータ学習部によって学習された前記第２確率分布のパラメータとに基づいて、前記第１確率分布と前記第２確率分布とを用いて表される規準を大きくするように、前記予測対象のソース音声に対応する前記ターゲット音声の各時刻の基本周波数を予測する
基本周波数パターン予測方法。
【請求項８】
請求項１〜請求項６の何れか１項に記載の基本周波数パターン予測装置の各部としてコンピュータを機能させるためのプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、基本周波数パターン予測装置、方法、及びプログラムに係り、特に、ソース音声から、ターゲット音声の基本周波数パターンを予測する基本周波数パターン予測装置、方法、及びプログラムに関する。
【背景技術】
【０００２】
他者とのコミュニケーションにおいて音声は利便性に優れた手段ではあるが、時として物理的制約により様々な障壁が必然的にもたらされる。例えば、発声器官の内、わずか一か所でも正常に動作しなくなると、深刻な発声障害を患い、音声コミュニケーションに支障をきたす。また、音声生成という物理的行為は、秘匿性の高い意思伝達には不向きであるし、周囲の騒音に脆弱である。これらの障壁を無くすためには、身体的制約を超えて発声器官を動作させて音声を生成したり、適切な発音動作を指定して音声を生成したり、聴取困難なほど微かな音声発声時の発声器官動作から通常音声を生成するなど、物理的・身体的制約を超えた音声生成機能の拡張が必要である。
【０００３】
例えば、喉頭癌などで喉頭を失った喉頭摘出者に対して、残存器官を用いた代替発声法により生成される自然性に乏しい音声を、より自然な音声へと変換する発声補助技術が提案されている（非特許文献１〜非特許文献３を参照）。この他にも、非可聴つぶやき音声を自然な音声に変換する技術も提案されており、秘匿性に優れた通話技術としての応用が期待されている。上述の技術はいずれも音声のスペクトル特徴量系列から自然音声の基本周波数(F0) パターンを予測する問題を扱っている点で共通しており、学習処理と変換処理で構成される。学習処理では、対象音声（前者であれば電気音声、後者であれば非可聴つぶやき音声）と通常音声の同一発話データを用いる。まず各離散時刻（以後，フレーム）において、前後数フレームから得られる対象音声のスペクトル特徴量と、通常音声の対数F0とその動的成分（時間微分または時間差分）を抽出し、スペクトル距離尺度に基づく動的時間伸縮によりこれらを対応付けた結合ベクトルを得る。これをパラレルデータと呼ぶ。各フレームのパラレルデータを用い、対象音声のスペクトル特徴量と通常音声の対数F0 の静的・動的成分の結合確率密度関数を混合正規分布モデル（Gaussian Mixture Model; GMM）で表現する。GMM のパラメータはExpectation-Maximization アルゴリズムにより学習することができる。変換処理では、学習されたGMM を用いて、系列内変動を考慮した最尤系列変換法により、対象音声のスペクトル特徴量系列から通常音声のF0 パターンへと変換することができる。
【先行技術文献】
【非特許文献】
【０００４】
【非特許文献１】Keigo Nakamura, Tomoki Toda, Hiroshi Saruwatari, Kiyohiro Shikano, "Speaking-aid systems using GMM-based voice conversion for electrolaryngeal speech," Speech Communication, vol. 54, no. 1, pp. 134-146, 2012.
【非特許文献２】Kou Tanaka, Tomoki Toda, Graham Neubig, Sakriani Sakti, Satoshi Nakamura, "A hybrid approach to electrolaryngeal speech enhancement based on noise reduction and statistical excitation generation," IEICE Transactions on Information and Systems, vol. E97-D, no. 6, pp. 1429-1437, Jun. 2014.
【非特許文献３】Kou Tanaka, Tomoki Toda, Graham Neubig, Sakriani Sakti, Satoshi Nakamura, "Direct F0 controlof an electrolarynx based on statistical excitation feature prediction and its evaluation through simulation," Proc. INTERSPEECH, pp. 31-35, Sep. 2014.
【非特許文献４】Hirokazu Kameoka, Jonathan Le Roux, Yasunori Ohishi, "A statistical model of speech F0 contours," ISCA Tutorial and Research Workshop on Statistical And Perceptual Audition (SAPA 2010), pp. 43-48, Sep. 2010.
【非特許文献５】Kota Yoshizato, Hirokazu Kameoka, Daisuke Saito, Shigeki Sagayama, "Hidden Markov convolutive mixture model for pitch contour analysis of speech," in Proc. The 13th Annual Conference of the International Speech Communication Association (Interspeech 2012), Sep. 2012.
【発明の概要】
【発明が解決しようとする課題】
【０００５】
従来技術では、学習処理や変換処理において音声のF0 パターンの物理的な生成過程を考慮したモデルが用いられていなかったため、物理的に人間が発声しえないような不自然なF0 パターンを生成することが起こりえた。この問題に対し、F0パターンの物理的な生成過程を考慮した予測を行うことで、より自然なF0 パターンを生成できる可能性がある。
【０００６】
F0パターンは声帯に張力を与える甲状軟骨の運動によって生み出されており、非特許文献４、５ではその制御機構の確率モデルに基づき、フレーズ・アクセント指令と呼ぶ甲状軟骨の運動に関係するパラメータを推定する技術が提案されている。この技術では，フレーズ・アクセント指令の時系列の生成プロセスを隠れマルコフモデル(HMM) により表現した点がポイントの一つであり，HMM のトポロジーの設計や遷移確率の学習を通して、指令列に関する言語学的ないし先験的な知識をパラメータ推定に組み込むことが可能である。
【０００７】
本発明は、上記事情を鑑みてなされたものであり、F0 パターンの物理的な生成過程の制約を考慮しながらスペクトル特徴量系列に対応する最適なF0パターンを推定することができる基本周波数パターン予測装置、方法、及びプログラムを提供することを目的とする。
【課題を解決するための手段】
【０００８】
上記の目的を達成するために本発明に係る基本周波数パターン予測装置は、学習サンプルのソース音声の時系列データとターゲット音声の時系列データとからなるパラレルデータを入力として、前記ソース音声の時系列データから抽出される各時刻のスペクトル特徴量ベクトルと、前記ターゲット音声の時系列データから抽出される、各時刻の基本周波数とに基づいて、前記ソース音声の各時刻のスペクトル特徴量ベクトルと、前記ターゲット音声の各時刻の基本周波数との間の関係をモデル化した第１確率分布のパラメータを学習する第１モデルパラメータ学習部と、前記ターゲット音声の各時刻の基本周波数に基づいて、基本周波数パターン生成過程をモデル化した第２確率分布のパラメータを学習する第２モデルパラメータ学習部と、予測対象のソース音声の時系列データを入力として、前記予測対象のソース音声の時系列データから抽出される各時刻のスペクトル特徴量ベクトルと、前記第１モデルパラメータ学習部によって学習された前記第１確率分布のパラメータと、前記第２モデルパラメータ学習部によって学習された前記第２確率分布のパラメータとに基づいて、前記第１確率分布と前記第２確率分布とを用いて表される規準を大きくするように、前記予測対象のソース音声に対応する前記ターゲット音声の各時刻の基本周波数を予測する基本周波数予測部と、を含んで構成されている。
【０００９】
本発明に係る基本周波数パターン予測方法は、第１モデルパラメータ学習部と、第２モデルパラメータ学習部と、基本周波数予測部とを含む基本周波数パターン予測装置における基本周波数パターン予測方法であって、前記第１モデルパラメータ学習部が、学習サンプルのソース音声の時系列データとターゲット音声の時系列データとからなるパラレルデータを入力として、前記ソース音声の時系列データから抽出される各時刻のスペクトル特徴量ベクトルと、前記ターゲット音声の時系列データから抽出される、各時刻の基本周波数とに基づいて、前記ソース音声の各時刻のスペクトル特徴量ベクトルと、前記ターゲット音声の各時刻の基本周波数との間の関係をモデル化した第１確率分布のパラメータを学習し、前記第２モデルパラメータ学習部が、前記ターゲット音声の各時刻の基本周波数に基づいて、基本周波数パターン生成過程をモデル化した第２確率分布のパラメータを学習し、前記基本周波数予測部が、予測対象のソース音声の時系列データを入力として、前記予測対象のソース音声の時系列データから抽出される各時刻のスペクトル特徴量ベクトルと、前記第１モデルパラメータ学習部によって学習された前記第１確率分布のパラメータと、前記第２モデルパラメータ学習部によって学習された前記第２確率分布のパラメータとに基づいて、前記第１確率分布と前記第２確率分布とを用いて表される規準を大きくするように、前記予測対象のソース音声に対応する前記ターゲット音声の各時刻の基本周波数を予測する。
【００１０】
本発明に係るプログラムは、上記の基本周波数パターン予測装置の各部としてコンピュータを機能させるためのプログラムである。
【発明の効果】
【００１１】
以上説明したように、本発明の基本周波数パターン予測装置、方法、及びプログラムによれば、ソース音声の各時刻のスペクトル特徴量ベクトルと、ターゲット音声の各時刻の基本周波数との間の関係をモデル化した第１確率分布のパラメータを学習し、基本周波数パターン生成過程をモデル化した第２確率分布のパラメータを学習し、予測対象のソース音声の時系列データから抽出される各時刻のスペクトル特徴量ベクトルから、第１確率分布と第２確率分布とを用いて表される規準を大きくするように、予測対象のソース音声に対応するターゲット音声の各時刻の基本周波数を予測することにより、F0 パターンの物理的な生成過程の制約を考慮しながらスペクトル特徴量系列に対応する最適なF0パターンを推定することができる、という効果が得られる。
【図面の簡単な説明】
【００１２】
【図１】ＨＭＭの状態遷移ネットワークの一例を説明するための図である。
【図２】ＨＭＭの状態遷移ネットワークの一例を説明するための図である。
【図３】ＨＭＭの状態遷移ネットワークの一例を説明するための図である。
【図４】状態の分割を説明するための図である。
【図５】本発明の第１の実施の形態に係る基本周波数パターン予測装置の構成を示す概略図である。
【図６】本発明の第１の実施の形態に係る基本周波数パターン予測装置の学習部の構成を示す概略図である。
【図７】本発明の第１の実施の形態に係る基本周波数パターン予測装置の変換処理部の構成を示す概略図である。
【図８】本発明の第１の実施の形態に係る基本周波数パターン予測装置における学習処理ルーチンの内容を示すフローチャートである。
【図９】本発明の第１の実施の形態に係る基本周波数パターン予測装置における基本周波数パターン予測処理ルーチンの内容を示すフローチャートである。
【図１０】本発明の第２の実施の形態に係る基本周波数パターン予測装置の変換処理部の構成を示す概略図である。
【図１１】本発明の第２の実施の形態に係る基本周波数パターン予測装置における学習処理ルーチンの内容を示すフローチャートである。
【図１２】本発明の第２の実施の形態に係る基本周波数パターン予測装置における基本周波数パターン予測処理ルーチンの内容を示すフローチャートである。
【図１３】実験データに用いた音声のＦ0パターンを示す図である。
【図１４】実験結果を示す図である。
【発明を実施するための形態】
【００１３】
以下、図面を参照して本発明の実施の形態を詳細に説明する。本発明で提案する技術は、音声の特徴量系列から基本周波数パターンを予測し、原音声の基本周波数パターンを予測した基本周波数パターンに置き換えることで音声の自然性を向上させることを目的とした音声処理技術である。
【００１４】
＜関連技術１：スペクトル特徴量系列からのF0 パターン予測方法＞
まず、スペクトル特徴量系列からのF0 パターン予測方法について説明する。
【００１５】
上記非特許文献１〜非特許文献３では、スペクトル特徴量系列からF0 パターンを予測する方法が提案されている。当該従来手法はスペクトル特徴量系列とF0 パターンの同時確率分布モデルのパラメータを学習する処理と学習した当該モデルを用いて所与のスペクトル特徴量系列からF0 パターンに変換する処理からなる。
【００１６】
＜学習処理＞
ソース音声（例えば電気音声）とターゲット音声（例えば自然音声）のパラレルデータが与えられているものとする。ソース音声のスペクトル特徴量ベクトルをc[k] とし，ターゲット音声の対数F0 とその動的成分（時間微分または時間差分）の結合ベクトル（F0 特徴量と呼ぶ。）をq[k] = (y[k];Δy[k])T とする。ここでk は離散時刻のインデックスである。音声特徴量c[k] としては例えば時刻k を中心とした前後数フレーム分のメルケプストラム（ベクトル）の系列を連結したベクトルに対し主成分分析により次元圧縮を行ったものを用いる。本手法ではc[k] とq[k] の同時確率分布を混合正規分布モデル（Gaussian Mixture Model; GMM）
【００１７】
【数１】
0006468519.tif 000002
【００１８】
でモデル化し、学習処理では所与のパラレルデータ｛c[k]; q[k]｝Kk=1 から当該GMM のパラメータ（各正規分布の重み，平均、分散共分散行列）を学習する。ただし、N(x;μ，Σ) は、xの確率密度関数が平均がμ、分散共分散行列がΣの正規分布で与えられることを意味する。
【００１９】
GMM のパラメータはExpectation-Maximization (EM) アルゴリズムにより推定することができる。学習したGMM パラメータをγとすると，条件付分布Ｐ(q[k]｜c[k]，γ) はスペクトル特徴量c[k] からF0 特徴量q[k] を予測するための分布と見ることができ、
【００２０】
【数２】
0006468519.tif 000003
【００２１】
のようにＰ(c[k]，q[k]｜γ) と同様にGMM で与えられる。ただし、
【００２２】
【数３】
0006468519.tif 000004
【００２３】
でありe(q｜c)mおよびD(q｜c)m は
【００２４】
【数４】
0006468519.tif 000005
【００２５】
で与えられる。
【００２６】
＜変換処理＞
変換処理では。所与のスペクトル特徴量系列
【００２７】
【数５】
0006468519.tif 000006
【００２８】
の下で，最尤のF0パターン
【００２９】
【数６】
0006468519.tif 000007
【００３０】
を以下の式（８）により求めることが目的である。
【００３１】
【数７】
0006468519.tif 000008
【００３２】
ただし、
【００３３】
【数８】
0006468519.tif 000009
【００３４】
であり、Ｗはｙとｑの関係を記述した変換行列（定数）である。ここで、P(ｑ｜ｃ，γ) は学習処理により学習したパラメータγで与えられるGMM であり、
【００３５】
【数９】
0006468519.tif 000010
【００３６】
で与えられる。ただし、ｍ= (ｍ1，...，ｍK) であり、ｍk は時刻k におけるGMMの成分インデックスを表す。ここで、P(q[k]｜c[k]，γ) が
【００３７】
【数１０】
0006468519.tif 000011
【００３８】
により近似できるとする。式(12)より、^mkは、データｃ[k]を生成したらしい確率が最も高い正規分布のインデックスを意味する。よって、式(9)よりP(ｑ｜ｃ，ｍ，＾γ）はすべてのkについて
【００３９】
【数１１】
0006468519.tif 000012
【００４０】
の積をとったもので与えられる。式(11)の近似によりP(ｑ｜ｃ，γ）は
【００４１】
【数１２】
0006468519.tif 000013
【００４２】
を連結したベクトルｅ(ｑ｜ｃ)を平均、
【００４３】
【数１３】
0006468519.tif 000014
【００４４】
を対角成分としたブロック対角行列Ｄ(q｜c) を分散共分散行列とした正規分布
【００４５】
【数１４】
0006468519.tif 000015
【００４６】
となる。これにｑ＝Ｗｙを代入し、ｙの分布となるように正規化すると、
【００４７】
【数１５】
0006468519.tif 000016
【００４８】
となるため、
【００４９】
【数１６】
0006468519.tif 000017
【００５０】
のようにｙの条件付き分布を得る。よって、式(8)の解は、
【００５１】
【数１７】
0006468519.tif 000018
【００５２】
となる。
【００５３】
＜関連技術２：F0パターン生成過程モデル＞
次に、F0パターン生成過程の確率モデルについて説明する。
【００５４】
音声のF0パターンの生成過程を記述したモデルに、藤崎の基本周波数（F0）パターン生成過程モデル（藤崎モデル）が知られている（非特許文献６）。
【００５５】
［非特許文献６］：H. Fujisaki, "In Vocal Physiology: Voice Production, Mechanisms and Functions," Raven Press, 1988.
【００５６】
藤崎モデルとは、甲状軟骨の運動によるF0 パターンの生成過程を説明した物理モデルである。藤崎モデルでは、甲状軟骨の二つの独立な運動（平行移動運動と回転運動）にそれぞれ伴う声帯の伸びの合計がF0の時間的変化をもたらすと解釈され、声帯の伸びとF0パターンの対数値y(t) が比例関係にあるという仮定に基づいてF0パターンがモデル化される。甲状軟骨の平行移動運動によって生じるF0パターンxp(t)をフレーズ成分、回転運動によって生じるF0 パターンxa(t) をアクセント成分と呼ぶ。藤崎モデルでは、音声のF0 パターンy(t) は、これらの成分に声帯の物理的制約によって決まるベースライン成分b を足し合わせたものとして、
【００５７】
【数１８】
0006468519.tif 000019
【００５８】
と表現される。これら二つの成分は二次の臨界制動系の出力と仮定され、
【００５９】
【数１９】
0006468519.tif 000020
【００６０】
と表される（＊ は時刻t に関する畳み込み演算）。ここでup(t) はフレーズ指令関数と呼ばれ、デルタ関数（フレーズ指令）の列からなり、ua(t) はアクセント指令関数と呼ばれ、矩形波（アクセント指令）の列からなる。これらの指令列には、発話の最初にはフレーズ指令が生起する、フレーズ指令は二連続で生起しない、異なる二つの指令は同時刻に生起しない、という制約条件がある。またαとβはそれぞれフレーズ制御機構、アクセント制御機構の固有角周波数であり、話者や発話内容によらず、おおよそα=3 rad/s、β=20 rad/s 程度であることが経験的に知られている。
【００６１】
＜関連技術３：F0パターン生成過程モデルパラメータ推定法＞
上述の藤崎モデルは以下のような確率モデルで記述することができる（非特許文献４、５、７参照）。
【００６２】
[非特許文献７]：石原達馬, 吉里幸太, 亀岡弘和, 齋藤大輔, 嵯峨山茂樹, \音声基本周波数の藤崎モデル指令列の統計的語彙モデル," 日本音響学会2013 年春季研究発表会講演論文集, 1-7-9, pp. 283-286, Mar. 2013.
【００６３】
まずフレーズ、アクセント指令関数のペアo[k] = (up[k], ua[k])T を出力するHMM を考える。ただし、k は離散時刻のインデックスを表す。状態出力分布は正規分布とし、各時刻の状態が与えられた下で
【００６４】
【数２０】
0006468519.tif 000021
【００６５】
により指令関数ペアo[k] が生成されるものとする。ここで｛sk｝Kk=1 はHMMの状態系列であり、平均ベクトルρ[k] はHMM の状態遷移の結果として定まる値である。具体的なHMM の構成の例を図１〜図３に示す。
【００６６】
図1に示すHMM の状態遷移ネットワークの例では、状態t = r0 においてμ(p)t [k] とμ(a)t はいずれも0 である。状態t = r0 からは状態p0 にのみ遷移することができ、状態t = p0 においてμ(p)t [k] は非負値A(p)[k] をとり、μ(a)t は0 となる。状態t = p0の次は状態r1 にのみ遷移することが許される。状態t = r0 同様、状態t = r1 においてμ(p)t [k] とμ(a)t はいずれも0 である。状態t = r1 からは状態a0，...，aN-1のいずれかにのみ遷移することができ、状態t = an においてμ(a)t は非負値A(a)nをとり、μ(p)t [k] は0 となる。状態t = anの次は状態r0 またはr1にのみ遷移することが許される。これよりμa[k] が矩形パルス列となることが保証される。
【００６７】
図2に示すHMMの状態遷移ネットワークの例では、状態t = r0 においてμ(p)t [k] とμ(a)t はいずれも0 である。状態t = r0 からは状態p0，...，pM-1のいずれかにのみ遷移することができ、状態t = pm においてμ(p)t は非負値A(p)mをとり、μ(a)t は0 となる。状態t = pmの次は状態r1にのみ遷移することが許される。状態t = r0 同様、状態t = r1 においてμ(p)t [k] とμ(a)t はいずれも0 である。状態t = r1 からは状態a0，...，aN-1のいずれかにのみ遷移することができ、状態t = an においてμ(a)t は非負値A(a)nをとり、μ(p)t は0 となる。状態t = anの次は状態r0 またはr1 にのみ遷移することが許される。これよりμa[k] が矩形パルス列となることが保証される。
【００６８】
図3に示すHMMの状態遷移ネットワークの例では、それぞれの終点と始点が連結された複数のLeft-to-Right 型HMM からなる。図１、２と同様、状態t = rlにおいてμ(p)t [k] とμ(a)t はいずれも0 である。また、状態t = pm においてμ(p)tは非負値A(p)m をとり、μ(a)t は0 となる。状態t = anにおいてμ(a)tは非負値A(a)nをとり、μ(p)tは0となる。
【００６９】
ρ[k]は、図1 の例では、以下の式（２４）で表わされる。
【００７０】
【数２１】
0006468519.tif 000022
【００７１】
また、図2, 3 の例では、ρ[k]は、以下の式（２５）で表わされる。
【００７２】
【数２２】
0006468519.tif 000023
【００７３】
いずれの例においても、図４のようにそれぞれの状態を同じ出力分布をもついくつかの小状態に分割し、Left-to-Right 型の状態遷移経路を制約することで同一状態に停留する時間長の確率をパラメータ化することができる。図4 は状態an を分割した例である。例えばこの図のように全てのm≠0 に対してan,n′からan,n′+1 への状態遷移確率を1に設定することで，an,0からan,n′への遷移確率が状態an がn′ステップだけ持続する確率に対応し、アクセント指令の持続長の確率を設定したり学習したりできるようになる。同様にpm とrl も小状態に分割することで、フレーズ指令の持続長と指令間の間隔の長さの分布をパラメータ化することが可能になる。以後、状態集合を
【００７４】
【数２３】
0006468519.tif 000024
【００７５】
と表記する。上記のHMM の構成は次のように書ける。
【００７６】
【数２４】
0006468519.tif 000025
【００７７】
状態系列s = ｛sk｝Kk=1 が与えられたとき、このHMM はフレーズ指令関数up[k] とアクセント指令関数ua[k] のペアを出力する。式(18) と式(20) で示した通り、up[k] とua[k] にそれぞれgp[k] とga[k] が畳み込まれてフレーズ成分xp[k] とアクセント成分xa[k] が出力される。これを式で表すと、
【００７８】
【数２５】
0006468519.tif 000026
【００７９】
と書ける（＊ は離散時刻k に関する畳み込み演算）。このとき，F0 パターンx[k] は
【００８０】
【数２６】
0006468519.tif 000027
【００８１】
と三種類の成分の重ね合わせで書ける。ただしb は時刻によらないベースライン成分である。
【００８２】
また、実音声においては、いつも信頼のできるF0 の値が観測できるとは限らない。藤崎モデルのパラメータ推定を行うにあたっては、信頼のおける観測区間のF0 値のみを考慮に入れて、そうでない区間は無視することが望ましい。例えば音声の無声区間においては通常声帯の振動に伴う周期的な粗密波は観測されないので、仮に自動ピッチ抽出によって音声の無声区間から何らかの値がF0 の推定値として得られたとしても、その値を声帯から発せられる信号のF0 の値と見なすのは適当ではない。そこで、提案モデルに観測F0値の時刻k における不確かさの程度ｖ2n[k] を導入する。具体的には、観測F0 値y[k] を、真のF0 値x[k] とノイズ成分
【００８３】
【数２７】
0006468519.tif 000028
【００８４】
との重ね合わせで
【００８５】
【数２８】
0006468519.tif 000029
【００８６】
と表現することで、信頼のおける区間かどうかに関わらず全ての観測区間を統一的に扱える。
【００８７】
xn[k] を周辺化することで、出力値系列o = ｛o[k]｝Kk=1 が与えられたときのy = ｛y[k]｝Kk=1 の確率密度関数
【００８８】
【数２９】
0006468519.tif 000030
【００８９】
が得られる。状態系列s = ｛sk｝Kk=1 と指令の振幅を表すパラメータ
【００９０】
【数３０】
0006468519.tif 000031
【００９１】
および遷移確率行列φ = (φi,j)I×I が与えられたとき、出力値系列o は
【００９２】
【数３１】
0006468519.tif 000032
【００９３】
に従って生成される。また、P(s｜φ) は状態遷移確率の積として
【００９４】
【数３２】
0006468519.tif 000033
と書ける。ただし、
【００９５】
【数３３】
0006468519.tif 000034
【００９６】
は初期状態がs1である確率をあらわす。式(30)、(32) および式(33) よりP(y,o,s｜θ,φ) は
【００９７】
【数３４】
0006468519.tif 000035
【００９８】
と書ける。これをo に関して周辺化すると
【００９９】
【数３５】
0006468519.tif 000036
【０１００】
が得られる。ただし、
【０１０１】
【数３６】
0006468519.tif 000037
0006468519.tif 000038
0006468519.tif 000039
【０１０２】
である。一方、s に関して周辺化すると
【０１０３】
【数３７】
0006468519.tif 000040
が得られる。ただし、Σs はあらゆる状態系列に関して和をとる操作を意味する。
【０１０４】
＜パラメータ推定アルゴリズム１＞
ｙとｏを完全データと見なすことで、式(35) を局所最大化するｓとθ をExpectation-Maximization アルゴリズムにより探索することができる。導出は省略するが、
【０１０５】
【数３８】
0006468519.tif 000041
【０１０６】
が大きくなるようにｓとθを更新するステップと、更新したｓとθを用いて
【０１０７】
【数３９】
0006468519.tif 000042
【０１０８】
とＲを
【０１０９】
【数４０】
0006468519.tif 000043
【０１１０】
により更新するステップを繰り返すことで式(35) を単調増加させることができる（詳細は、上記非特許文献４参照）。
【０１１１】
具体的には、以下の初期設定、Ｅステップ、及びＭステップが実行される。
【０１１２】
（初期設定）
ｓとθを初期設定する。
【０１１３】
（E ステップ）
フレーズ成分、アクセント成分、ベースライン成分の条件付き期待値
【０１１４】
【数４１】
0006468519.tif 000044
【０１１５】
と条件付き分散共分散行列Ｒを
【０１１６】
【数４２】
0006468519.tif 000045
により更新する。ただし、
【０１１７】
【数４３】
0006468519.tif 000046
0006468519.tif 000047
【０１１８】
である。また、Ｒにおける
【０１１９】
【数４４】
0006468519.tif 000048
【０１２０】
に対応するブロック対角成分を
【０１２１】
【数４５】
0006468519.tif 000049
【０１２２】
とする。
【０１２３】
すなわち、
【０１２４】
【数４６】
0006468519.tif 000050
【０１２５】
である（＊ は以後用いないブロック成分である）。
【０１２６】
（Ｍステップ）
Q(ｓ,θ) が最大となる状態系列ｓ= (ｓ1,...,ｓK) を探索する。γp とγaは対角行
【０１２７】
列なので、
【０１２８】
【数４７】
0006468519.tif 000051
0006468519.tif 000052
【０１２９】
はいずれも
【０１３０】
【数４８】
0006468519.tif 000053
【０１３１】
のようにｋごとの項の和の形で書ける。従って、Ｑ(ｓ，θ) はｓに依らない項を除けば
【０１３２】
【数４９】
0006468519.tif 000054
【０１３３】
と書ける。従って、Ｑ(ｓ,θ) を最大にする状態系列ｓ= (ｓ1,...,ｓK)はViterbiアルゴリズムにより求めることができる（詳細は上記非特許文献４参照）。ただし、[・]k,k は行列のk 行k 列の要素、[・]kはベクトルの第k 要素を表す。
【０１３４】
続いて、Ｑ(ｓ,θ)を最大にするようにθを更新する。Ｑ(ｓ,θ)を最大にするθは、Ｑ(ｓ,θ)の各変数に関する偏微分を0 と置くことにより得られる（（詳細は上記非特許文献４参照）。
【０１３５】
また、推定された状態系列ｓから、状態遷移確率φが求められる。
【０１３６】
［第１の実施の形態］
＜本発明の実施の形態の概要＞
本発明の実施の形態の技術は、上述した関連技術１のスペクトル特徴量系列からのF0 パターン予測方法と同様、学習処理と変換処理からなるが、式(8) の代わりに、上述した関連技術１のスペクトル特徴量系列からのF0 パターン予測方法の確率分布と、上述した関連技術３のF0 パターン生成過程モデルパラメータ推定法の確率分布の積を規準とすることにより、上述したF0パターン生成過程モデルにできるだけ即したF0パターンをスペクトル特徴量から統計的に予測することを可能にする技術である。
【０１３７】
学習処理ではパラレルデータの学習サンプル｛c[k]，q[k]｝Kk=1 が与えられた下で
【０１３８】
【数５０】
0006468519.tif 000055
【０１３９】
ができるだけ大きくなるようにγを学習する。また、学習サンプルの基本周波数パターン｛y[k]｝Kk=1 が与えられた下でＰ(ｙ，ｓ｜θ、φ)ができるだけ大きくなるようにθ とφを学習する。
【０１４０】
一方、変換処理では入力音声の｛c[k]｝Kk=1 が与えられた下でＰ(q｜c，γ)Ｐ(ｙ，ｓ｜θ，φ)またはこれらを近似する確率密度関数ができるだけ大きくなるようにｙを求める。
【０１４１】
＜システム構成＞
次に、ソース音声のスペクトル特徴量系列から、ターゲット音声の基本周波数パターンを予測する基本周波数パターン予測装置に、本発明を適用した場合を例にして、本発明の実施の形態を説明する。
【０１４２】
図５に示すように、本発明の第１の実施の形態に係る基本周波数パターン予測装置は、ＣＰＵと、ＲＡＭと、後述する学習処理ルーチン、及び基本周波数パターン予測処理ルーチンを実行するためのプログラムを記憶したＲＯＭとを備えたコンピュータで構成され、機能的には次に示すように構成されている。
【０１４３】
図５に示すように、基本周波数パターン予測装置１００は、入力部１０と、演算部２０と、出力部９０とを備えている。
【０１４４】
入力部１０は、学習サンプルのソース音声（例えば電気音声）の時系列データとターゲット音声（例えば自然音声）の時系列データとからなるパラレルデータを受け付ける。また、入力部１０は、予測対象のソース音声の時系列データを受け付ける。
【０１４５】
演算部２０は、学習部３０と、パラメータ記憶部４０と、変換処理部５０とを備えている。
【０１４６】
図６に示すように、学習部３０は、特徴量抽出部３２と、基本周波数系列抽出部３４と、第１モデルパラメータ学習部３６と、第２モデルパラメータ学習部３８とを備えている。
【０１４７】
特徴量抽出部３２は、入力部１０によって受け付けた学習サンプルのソース音声の時系列データから、ソース音声のスペクトグラム特徴量ベクトルc[k]を抽出する。ここでk は離散時刻のインデックスである。例えば、非特許文献１〜３と同様に時刻k を中心とした前後数フレーム分のメルケプストラム（ベクトル）の系列を連結したベクトルに対し主成分分析により次元圧縮を行ったものをc[k] として用いる。
【０１４８】
基本周波数系列抽出部３４は、入力部１０によって受け付けた学習サンプルのターゲット音声の時系列データから、ターゲット音声の各時刻ｋにおける基本周波数y[k]を抽出し、y = (y[1]，...， y[K])Tとする。
【０１４９】
この基本周波数の抽出処理は、周知技術により実現でき、例えば、非特許文献８（H. Kameoka, "Statistical speech spectrum model incorporating all-pole vocal tract model and F0 contour generating process model," in Tech. Rep. IEICE, 2010, in Japanese.）に記載の手法を利用して、８ｍｓごとに基本周波数を抽出する。
【０１５０】
また、y とその動的成分（時間微分または時間差分）の結合ベクトル（F0 特徴量と呼ぶ。）をq[k] = (y[k]，Δy[k])T とする。
【０１５１】
以上より、｛c[k]，q[k]｝Kk=1 というデータが得られる。
【０１５２】
第１モデルパラメータ学習部３６は、特徴量抽出部３２によって抽出された各時刻ｋのスペクトル特徴量ベクトルc[k]と、基本周波数系列抽出部３４によって抽出された各時刻ｋの基本周波数の結合ベクトルｑ［ｋ］とに基づいて、ソース音声の各時刻のスペクトル特徴量ベクトルc[k]とターゲット音声の各時刻ｋの基本周波数の結合ベクトルｑ［ｋ］との同時確率分布を表す混合正規分布である第１確率分布のパラメータを学習する。
【０１５３】
具体的には、第１モデルパラメータ学習部３６は、上述したスペクトル特徴量系列からのＦ0パターン予測方法の学習処理と同様に、式(1) のGMM のパラメータγを学習する。学習したGMM パラメータを^γとする。
【０１５４】
第２モデルパラメータ学習部３８は、基本周波数系列抽出部３４によって抽出された各時刻ｋの基本周波数y[k]に基づいて、各時刻ｋの基本周波数y[k]と、隠れマルコフモデルの各時刻の状態からなる状態系列ｓとの組み合わせの確率分布である第２確率分布のパラメータを学習する。
【０１５５】
具体的には、第２モデルパラメータ学習部３８は、上述した関連技術３のF0 パターン生成過程モデルパラメータ推定法のパラメータ推定アルゴリズム１に従って、F0パターン生成過程モデルのパラメータθ、φを学習する。
【０１５６】
もし学習サンプルのフレーズ指令系列とアクセント指令系列のペアのデータo = ｛ok｝Kk=1 が入手できるのであれば，o からθ、φを学習しても良い（HMM の通常の学習に相当）。学習したF0パターン生成過程モデルのパラメータを^θ、^φとする。
【０１５７】
変換処理部５０は、予測対象のソース音声の時系列データを入力として、ソース音声の時系列データから抽出される各時刻のスペクトル特徴量ベクトルと、第１モデルパラメータ学習部３６によって学習された第１確率分布のパラメータγと、第２モデルパラメータ学習部３８によって学習された第２確率分布のパラメータθ、φとに基づいて、第１確率分布と第２確率分布との積を用いて表される規準を大きくするように、各時刻の基本周波数ｙと、各時刻の状態からなる状態系列ｓを推定することにより、予測対象のソース音声に対応するターゲット音声の各時刻の基本周波数ｙを予測する。
【０１５８】
ここで、予測対象のソース音声に対応するターゲット音声の各時刻の基本周波数ｙを予測する原理について説明する。
【０１５９】
学習処理で学習したパラメータ^γ、＾θ、^φと、予測対象のソース音声の特徴量系列c = ｛c[k]｝Kk=1 を用いて、式(9)と式(35) の積
【０１６０】
【数５１】
0006468519.tif 000056
ができるだけ大きくなるようにｙ、ｓを推定する。ただし、ｑ = Ｗｙであり、
【０１６１】
【数５２】
0006468519.tif 000057
【０１６２】
である。ωはF0パターンの予測においてF0パターン生成過程のモデルをどれだけ考慮に入れるかを意味した非負の定数である。
【０１６３】
以下に、
【０１６４】
【数５３】
0006468519.tif 000058
【０１６５】
を大きくするためのアルゴリズムについて述べる。上述した関連技術１のスペクトル特徴量系列からのF0 パターン予測方法と同様、
【０１６６】
【数５４】
0006468519.tif 000059
【０１６７】
と近似することで以下の反復処理によりｙ、ｓを推定することができる（ステップ１と２の実行順序は任意）。
【０１６８】
（ステップ１）上述した関連技術１のスペクトル特徴量系列からのF0 パターン予測方法の変換処理によりｙを初期設定する。
【０１６９】
（ステップ２）cを用いて^ｍを式(87) により求める。
【０１７０】
（ステップ３）ｙを固定し、上述した関連技術３のF0 パターン生成過程モデルパラメータ推定法のパラメータ推定アルゴリズム１によりｓを推定する。
【０１７１】
（ステップ４）ｓと^ｍを固定して以下の式によりｙを更新し、ステップ３に戻る。
【０１７２】
【数５５】
0006468519.tif 000060
【０１７３】
（ＧγＧT）-1は大きなサイズの行列の逆行列であるが、以下に示すやり方で効率的に計算することができる。ＧγＧTは
【０１７４】
【数５６】
0006468519.tif 000061
【０１７５】
であること、Ｇ-1p とＧ-1a が
【０１７６】
【数５７】
0006468519.tif 000062
【０１７７】
のような下三角帯行列で近似できることより、（ＧγＧT）-1は
【０１７８】
【数５８】
0006468519.tif 000063
【０１７９】
と書け、さらにWoodbury の公式
【０１８０】
【数５９】
0006468519.tif 000064
【０１８１】
を式(95) 右辺に適用することで（ＧγＧT）-1は
【０１８２】
【数６０】
0006468519.tif 000065
【０１８３】
と書ける。さらにWoodbury の公式より
【０１８４】
【数６１】
0006468519.tif 000066
【０１８５】
は
【０１８６】
【数６２】
0006468519.tif 000067
【０１８７】
と書ける。式(93)、(94) より
【０１８８】
【数６３】
0006468519.tif 000068
【０１８９】
はいずれも帯行列になるので
【０１９０】
【数６４】
0006468519.tif 000069
【０１９１】
の形の計算はCholesky 分解により効率的に計算することができる。ただし、a は任意のベクトル，A は任意の行列である。
【０１９２】
以上説明した原理を実現するために、本実施の形態では、図７に示すように、変換処理部５０は、特徴量抽出部５２と、基本周波数系列予測部５４と、正規分布系列予測部５６と、状態系列推定部５８と、基本周波数系列更新部６０と、収束判定部６２とを備えている。
【０１９３】
特徴量抽出部５２は、入力部１０によって受け付けた予測対象のソース音声の時系列データから、特徴量抽出部３２と同様に、ソース音声の各時刻ｋのスペクトグラム特徴量ベクトルc[k]を抽出する。
【０１９４】
基本周波数系列予測部５４は、第１モデルパラメータ学習部３６によって学習された第１確率分布のパラメータγと、特徴量抽出部５２によって抽出された各時刻ｋのスペクトル特徴量ベクトルc[k]とに基づいて、上述したＦ0パターン予測方法の変換処理と同様に、上記式（１６）に従って、各時刻ｋの基本周波数ｙ［ｋ］を推定することにより、各時刻ｋの基本周波数ｙ［ｋ］を初期設定する。
【０１９５】
正規分布系列予測部５６は、第１モデルパラメータ学習部３６によって学習された第１確率分布のパラメータγと、特徴量抽出部５２によって抽出された各時刻ｋのスペクトル特徴量ベクトルc[k]とに基づいて、上記式（８７）に従って、各時刻ｋのスペクトル特徴量ベクトルc[k]を生成したらしい確率が最も高い正規分布のインデックス＾ｍkを推定する。
【０１９６】
状態系列推定部５８は、基本周波数系列予測部５４によって初期設定された、または状態系列推定部５８によって前回更新された各時刻ｋの基本周波数ｙ［ｋ］を固定して、上述した関連技術３のF0 パターン生成過程モデルパラメータ推定法のパラメータ推定アルゴリズム１と同様に、上記式（３５）を局所最大化する状態系列ｓと各時刻ｋにおける状態に応じたフレーズ指令の振幅及び各アクセント指令の振幅を表すパラメータθとを、ＥＭアルゴリズムにより探索することにより、状態系列ｓを推定する。
【０１９７】
基本周波数系列更新部６０は、状態系列推定部５８によって推定された状態系列ｓと、正規分布系列予測部５６によって推定された各時刻の正規分布のインデックス＾ｍkとに基づいて、上記式（８８）に従って、各時刻ｋの基本周波数ｙ［ｋ］を更新する。
【０１９８】
収束判定部６２は、予め定められた収束判定条件を満たすまで、状態系列推定部５８及び基本周波数系列更新部６０による各処理を繰り返させる。収束判定条件としては、例えば、予め定められた繰り返し回数に到達することである。
【０１９９】
収束判定条件を満たしたときに、最終的に得られた各時刻ｋの基本周波数ｙ［ｋ］を、予測対象のソース音声に対応するターゲット音声の各時刻の基本周波数の予測結果として、出力部９０により出力する。
【０２００】
＜基本周波数パターン予測装置の作用＞
次に、本実施の形態に係る基本周波数パターン予測装置１００の作用について説明する。まず、学習サンプルのソース音声の時系列データ及びターゲット音声の時系列データからなるパラレルデータが、基本周波数パターン予測装置１００に入力されると、基本周波数パターン予測装置１００において、図８に示す学習処理ルーチンが実行される。
【０２０１】
まず、ステップＳ１０１において、入力されたソース音声の時系列データを読み込み、
【０２０２】
各時刻ｋのスペクトル特徴量ベクトルｃ［ｋ］を抽出する。ステップＳ１０２において、入力されたターゲット音声の時系列データを読み込み、ターゲット音声の各時刻ｋにおける基本周波数y[k]を抽出し、また、基本周波数y[k]とその動的成分の結合ベクトルq[k]を抽出する。
【０２０３】
そして、ステップＳ１０３において、上記ステップＳ１０１で抽出された各時刻ｋのスペクトル特徴量ベクトルc[k]と、上記ステップＳ１０２で抽出された各時刻ｋの基本周波数の結合ベクトルｑ［ｋ］とに基づいて、上記式(1) のGMM のパラメータγを学習する。
【０２０４】
ステップＳ１０４では、状態系列ｓと、各時刻における状態に応じたフレーズ指令の振幅及び各アクセント指令の振幅を表すパラメータθとを初期設定する。
【０２０５】
そして、ステップＳ１０５において、上記式（４９）、式（５０）に従って、フレーズ成分、アクセント成分、ベースライン成分の条件付き期待値￣ｘと、条件付き分散共分散行列Ｒとを更新する。
【０２０６】
次のステップＳ１０６では、上記ステップＳ１０４で初期設定された、又は後述するステップＳ１０７で前回更新されたパラメータθと、上記ステップＳ１０５で更新されたフレーズ成分、アクセント成分、ベースライン成分の条件付き期待値￣ｘと、条件付き分散共分散行列Ｒとに基づいて、上記式（６３）式を用いて、Ｑ(ｓ,θ) を最大にする状態系列ｓ= (ｓ1,...,ｓK)をViterbiアルゴリズムにより求めて、状態系列ｓを更新する。
【０２０７】
ステップＳ１０７では、上記ステップＳ１０６で更新された状態系列ｓと、上記ステップＳ１０５で更新されたフレーズ成分、アクセント成分、ベースライン成分の条件付き期待値￣ｘと、条件付き分散共分散行列Ｒとに基づいて、Ｑ(ｓ,θ)の各変数に関する偏微分を0 と置くことにより、Ｑ(ｓ,θ) を最大にするパラメータθを求めて、各時刻における状態に応じたフレーズ指令の振幅及び各アクセント指令の振幅を表すパラメータθとを更新する。
【０２０８】
ステップＳ１０８において、予め定められた収束判定条件を満たしたか否かを判定し、収束判定条件を満たしていない場合には、上記ステップＳ１０５へ戻る。一方、収束判定条件を満たした場合には、ステップＳ１０９において、上記ステップＳ１０３で学習されたパラメータγ、上記ステップＳ１０７で最終的に得られたパラメータθとを、パラメータ記憶部４０に格納する。また、上記ステップＳ１０６で最終的に得られた状態系列ｓから、状態遷移確率φを求め、パラメータ記憶部４０に格納する。
【０２０９】
次に、予測対象のソース音声の時系列データが、基本周波数パターン予測装置１００に入力されると、基本周波数パターン予測装置１００において、図９に示す基本周波数パターン予測処理ルーチンが実行される。
【０２１０】
まず、ステップＳ１２１において、入力された予測対象のソース音声の時系列データを読み込み、各時刻ｋのスペクトル特徴量ベクトルｃ［ｋ］を抽出する。ステップＳ１２２において、パラメータ記憶部４０に記憶されたパラメータγと、上記ステップＳ１２１で抽出された各時刻のスペクトル特徴量ベクトルc[k]とに基づいて、上記式（１６）に従って、各時刻ｋの基本周波数ｙ［ｋ］を推定することにより、各時刻ｋの基本周波数ｙ［ｋ］を初期設定する。
【０２１１】
そして、ステップＳ１２３では、パラメータ記憶部４０に記憶されたパラメータγと、上記ステップＳ１２１で抽出された各時刻ｋのスペクトル特徴量ベクトルc[k]とに基づいて、上記式（８７）に従って、各時刻ｋのスペクトル特徴量ベクトルc[k]を生成したらしい確率が最も高い正規分布のインデックス＾ｍkを推定する。
【０２１２】
ステップＳ１２４では、上記ステップＳ１２２で初期設定された、または後述するステップＳ１２５で前回更新された各時刻ｋの基本周波数ｙ［ｋ］を固定して、上記式（３５）を局所最大化する状態系列ｓと各時刻ｋにおける状態に応じたフレーズ指令の振幅及び各アクセント指令の振幅を表すパラメータθとを、ＥＭアルゴリズムにより探索することにより、状態系列ｓを推定する。
【０２１３】
ステップＳ１２５では、上記ステップＳ１２４で推定された状態系列ｓと、上記ステップＳ１２３で推定された各時刻の正規分布のインデックス＾ｍkとに基づいて、上記式（８８）に従って、各時刻ｋの基本周波数ｙ［ｋ］を更新する。
【０２１４】
ステップＳ１２６では、予め定められた収束判定条件を満たしたか否かを判定し、収束判定条件を満たしていない場合には、ステップＳ１２４へ戻る。一方、収束判定条件を満たした場合には、ステップＳ１２７において、上記ステップＳ１２５で最終的に得られた各時刻ｋの基本周波数ｙ［ｋ］を、予測対象のソース音声に対応するターゲット音声の各時刻の基本周波数の予測結果として、出力部９０により出力し、基本周波数パターン予測処理ルーチンを終了する。
【０２１５】
以上説明したように、第１の実施の形態に係る基本周波数パターン予測装置によれば、ソース音声の各時刻のスペクトル特徴量ベクトルと、ターゲット音声の各時刻の基本周波数との間の関係をモデル化したＧＭＭである第１確率分布Ｐ(q[k]，c[k]｜γ)のパラメータγを学習し、基本周波数パターン生成過程をモデル化した第２確率分布Ｐ（ｙ，ｓ｜θ，φ）のパラメータθ、φを学習し、予測対象のソース音声の時系列データから抽出される各時刻のスペクトル特徴量ベクトルから、第１確率分布Ｐ(q[k]，c[k]｜γ)と第２確率分布Ｐ（ｙ，ｓ｜θ，φ）との積を用いて表される規準を大きくするように、予測対象のソース音声に対応するターゲット音声の各時刻の基本周波数を予測することにより、F0 パターンの物理的な生成過程の制約を考慮しながらスペクトル特徴量系列に対応する最適なF0パターンを推定することができる。
【０２１６】
［第２の実施の形態］
次に、本発明の第２の実施の形態に係る基本周波数パターン予測装置について説明する。なお、第１の実施の形態と同様の構成となる部分については同一符号を付して説明を省略する。
【０２１７】
第２の実施の形態では、第２の確率分布及びパラメータを推定する方法と、各時刻の基本周波数の予測方法とが第１の実施の形態と異なっている。
【０２１８】
第２の実施の形態に係る基本周波数パターン予測装置の学習部３０の第２モデルパラメータ学習部３８によるパラメータを学習する原理について説明する。
【０２１９】
まず、関連技術３のF0パターン生成過程モデルパラメータ推定法のパラメータ推定アルゴリズムについて説明する。
【０２２０】
＜パラメータ推定アルゴリズム２＞
観測F0系列ｙが与えられたもとで、モデルパラメータθとｏの事後確率Ｐ(ｏ，θ｜ｙ)の局所最適解を求める反復アルゴリズムを以下に示す。状態系列ｓを隠れ変数とし、事後確率Ｐ(ｏ，θ｜ｙ) が、
【０２２１】
【数６５】
0006468519.tif 000070
【０２２２】
をｓについて周辺化することで得られる点に注意すると、Ｑ関数Ｑ(ｏ，θ，ｏ´，θ´)は
【０２２３】
【数６６】
0006468519.tif 000071
と置ける。ただし、
【０２２４】
【数６７】
0006468519.tif 000072
【０２２５】
は定数項を除いて等しいことを表す。また、ｇb[k] = δ[k] (クロネッカーのデルタ) である。よって、Ｐ(ｓk = ｔ｜ｙ，ｏ´，θ´)をForward-Backward アルゴリズムにより計算するステップ、ｏとθについてＱ(ｏ，θ，ｏ´，θ´)を増加させるステップを繰り返すことで、Ｐ(ｏ，θ｜ｙ) が局所最大となる解を得ることができる。ｏはフレーズ・アクセント指令系列のペアであるため、Ｑ(ｏ，θ，ｏ´，θ´)を増加させるステップにおいては、ｏの非負制約を考慮する必要がある。ｏの非負制約を満たしながらＱ(ｏ，θ，ｏ´，θ´)を増加させるような更新則は以下により導くことができる。まず、Ｑ(ｏ，θ，ｏ´，θ´)の下界はJensen の不等式より
【０２２６】
【数６８】
0006468519.tif 000073
【０２２７】
のように設計することができる。また、i，k，l は、
【０２２８】
【数６９】
0006468519.tif 000074
【０２２９】
を満たす任意の変数である。従ってQ関数の下界は、
【０２３０】
【数７０】
0006468519.tif 000075
0006468519.tif 000076
【０２３１】
と表される。この下界関数をλi,k,l≧0 に関して最大化するステップとｏ に関して最大化するステップを交互に繰り返せばＱ(ｏ，θ，ｏ´，θ´)を増加させることができる。いずれのステップの更新則も解析的に求めることができ、それぞれ
【０２３２】
【数７１】
0006468519.tif 000077
【０２３３】
で表される。以上の反復が収束したあと、続けてθを更新する。更新式は、図１の場合、
【０２３４】
【数７２】
0006468519.tif 000078
【０２３５】
図２、３の場合、
【０２３６】
【数７３】
0006468519.tif 000079
【０２３７】
である。これらの更新値をｏ´とθ´に代入したのちに、Ｐ(ｓk = ｔ｜ｙ，ｏ´，θ´)の更新を再度行い、以後同様の処理を繰り返すことで事後確率Ｐ(ｏ，θ｜ｙ)を単調増加させることができる。
【０２３８】
以上の反復アルゴリズムが収束した後、上述したパラメータ推定アルゴリズム１のViterbi アルゴリズムにより求まる最適なｓを状態系列の推定結果とする。
【０２３９】
また、推定された状態系列ｓから、状態遷移確率φが求められる。
【０２４０】
以上説明した原理に従って、第２モデルパラメータ学習部３８は、基本周波数系列抽出部３４によって抽出された各時刻ｋの基本周波数y[k]に基づいて、各時刻ｋの基本周波数y[k]と、各時刻ｋにおける甲状軟骨の平行移動運動によって生じる基本周波数パターンを表すフレーズ指令ｕp［ｋ］及び甲状軟骨の回転運動によって生じる基本周波数パターンを表すアクセント指令ｕa［ｋ］のペアからなる指令関数ｏ［ｋ］との組み合わせの確率分布である第２確率分布のパラメータθ、φを学習する。
【０２４１】
第２の実施の形態における変換処理部５０は、予測対象のソース音声の時系列データを入力として、ソース音声の時系列データから抽出される各時刻ｋのスペクトル特徴量ベクトルｃ［ｋ］と、第１モデルパラメータ学習部３６によって学習された第１確率分布のパラメータγと、第２モデルパラメータ学習部３８によって学習された第２確率分布のパラメータθ、φとに基づいて、第１確率分布と第２確率分布とを用いて表される規準を大きくするように、各時刻ｋの基本周波数ｙ［ｋ］と、各時刻ｋのフレーズ指令及びアクセント指令のペアからなる指令関数ｏ［ｋ］とを推定することにより、予測対象のソース音声に対応するターゲット音声の各時刻ｋの基本周波数ｙ［ｋ］を予測する。
【０２４２】
ここで、予測対象のソース音声に対応するターゲット音声の各時刻ｋの基本周波数ｙ［ｋ］を予測する原理について説明する。
【０２４３】
＜変換処理＞
学習処理で学習したパラメータ＾γ、＾θ、＾ψと、予測対象のソース音声の特徴量系列ｃ＝｛c[k]｝Kk=1 を用いて、式(9) と式(45) の積
【０２４４】
【数７４】
0006468519.tif 000080
【０２４５】
ができるだけ大きくなるようにｙ、ｏを推定する。ただし、ｑ = Ｗｙであり、
【０２４６】
【数７５】
0006468519.tif 000081
である。
【０２４７】
次に、
【０２４８】
【数７６】
0006468519.tif 000082
【０２４９】
を大きくするためのアルゴリズムについて述べる。上述した関連技術１のスペクトル特徴量系列からのF0 パターン予測方法と同様に、
【０２５０】
【数７７】
0006468519.tif 000083
【０２５１】
と近似することで以下の反復処理によりｙ、ｏを推定することができる（ステップ1 と2 の実行順序は任意）。
【０２５２】
（ステップ１）上述した関連技術１のスペクトル特徴量系列からのF0 パターン予測方法の変換処理によりｙを初期設定する。
【０２５３】
（ステップ２）cを用いて^ｍを式(１００) により求める。
【０２５４】
（ステップ３）ｙを固定し、上述した関連技術３のF0 パターン生成過程モデルパラメータ推定法のパラメータ推定アルゴリズム２によりｏを推定する。
【０２５５】
（ステップ４）ｏと＾ｍを固定して以下の式によりｙを更新し、ステップ３に戻る。
【０２５６】
【数７８】
0006468519.tif 000084
【０２５７】
以上説明した原理を実現するために、第２の実施の形態では、図１０に示すように、変換処理部５０は、特徴量抽出部５２と、基本周波数系列予測部５４と、正規分布系列予測部２５６と、指令系列推定部２５８と、基本周波数系列更新部２６０と、収束判定部６２とを備えている。
【０２５８】
正規分布系列予測部２５６は、第１モデルパラメータ学習部３６によって学習された第１確率分布のパラメータγと、特徴量抽出部５２によって抽出された各時刻のスペクトル特徴量ベクトルc[k]とに基づいて、上記式（１００）に従って、各時刻ｋのスペクトル特徴量ベクトルc[k]を生成したらしい確率が最も高い正規分布のインデックス＾ｍkを推定する。
【０２５９】
指令系列推定部２５８は、基本周波数系列予測部５４によって初期設定された、または状態系列推定部５８によって前回更新された各時刻ｋの基本周波数ｙ［ｋ］を固定して、上述した関連技術３のF0 パターン生成過程モデルパラメータ推定法のパラメータ推定アルゴリズム２と同様に、事後確率Ｐ（ｏ，θ｜ｙ）を局所最大化する指令系列ｏを推定する。
【０２６０】
基本周波数系列更新部２６０は、指令系列推定部２５８によって推定された指令系列ｏと、正規分布系列予測部５６によって推定された各時刻の正規分布のインデックス＾ｍkとに基づいて、上記式（１０１）に従って、各時刻ｋの基本周波数ｙ［ｋ］を更新する。
【０２６１】
収束判定部６２は、予め定められた収束判定条件を満たすまで、指令系列推定部２５８及び基本周波数系列更新部２６０による各処理を繰り返させる。収束判定条件としては、例えば、予め定められた繰り返し回数に到達することである。
【０２６２】
収束判定条件を満たしたときに、最終的に得られた各時刻ｋの基本周波数ｙ［ｋ］を、予測対象のソース音声に対応するターゲット音声の各時刻の基本周波数の予測結果として、出力部９０により出力する。
【０２６３】
＜基本周波数パターン予測装置の作用＞
【０２６４】
次に、第２の実施の形態に係る基本周波数パターン予測装置１００の作用について説明する。なお、第１の実施の形態と同様の処理については、同一符号を付して詳細な説明を省略する。
【０２６５】
まず、学習サンプルのソース音声の時系列データ及びターゲット音声の時系列データからなるパラレルデータが、基本周波数パターン予測装置１００に入力されると、基本周波数パターン予測装置１００において、図１１に示す学習処理ルーチンが実行される。
【０２６６】
まず、ステップＳ１０１において、入力されたソース音声の時系列データを読み込み、各時刻ｋのスペクトル特徴量ベクトルｃ［ｋ］を抽出する。ステップＳ１０２において、入力されたターゲット音声の時系列データを読み込み、ターゲット音声の各時刻ｋにおける基本周波数y[k]を抽出し、また、各時刻ｋの基本周波数y[k]とその動的成分の結合ベクトルq[k]を抽出する。
【０２６７】
そして、ステップＳ１０３において、上記ステップＳ１０１で抽出された各時刻ｋのスペクトル特徴量ベクトルc[k]と、上記ステップＳ１０２で抽出された各時刻ｋの基本周波数の結合ベクトルｑ［ｋ］とに基づいて、上記式(1) のGMM のパラメータγを学習する。
【０２６８】
ステップＳ２００では、指令系列ｏと、各時刻における状態に応じたフレーズ指令の振幅及び各アクセント指令の振幅を表すパラメータθとを初期設定する。また、ターゲット音声の時系列データに基づいて、有声区間、無声区間を特定し、各時刻ｋの基本周波数の不確かさの程度vn2 [ｋ]を推定する。
【０２６９】
そして、ステップＳ２０１において、上記ステップＳ２００で設定された指令系列oの初期値、または後述するステップＳ２０３で前回更新された指令系列oに基づいて、（k,t）の全ての組み合わせについて、事後確率Ｐ(sk=t|y,o′,θ′)を更新する。
【０２７０】
ステップＳ２０２では、上記ステップＳ２００で設定された指令系列oの初期値、または後述するステップＳ２０３で前回更新された指令系列oに基づいて、（ｋ、ｌ）の全ての組み合わせについて、上記の式（７１）に従って、補助変数λp,k,l、λa,k,l、λb,k,lを算出して更新する。
【０２７１】
次のステップＳ２０３では、上記ステップＳ１０２で抽出されたされた基本周波数系列ｙと、上記ステップＳ２００で算出された各時刻ｋの不確かさの程度vn2 [ｋ]と、上記ステップＳ２０１で更新された事後確率Ｐ(sk=t|y,o′,θ′)と、上記ステップＳ２０２で更新された補助変数λp,k,l、λa,k,l、λb,k,lとに基づいて、上記式（７２）に従って、非負値である各時刻ｌのフレーズ指令ｕp［ｌ］及びアクセント指令ｕa［ｌ］からなる指令系列oとベース成分ｕbとを更新する。
【０２７２】
次のステップＳ２０４では、収束条件として、繰り返し回数ｓが、Ｓに到達したか否かを判定し、繰り返し回数ｓがＳに到達していない場合には、収束条件を満足していないと判断して、上記ステップＳ２０２へ戻る。一方、繰り返し回数ｓがＳに到達した場合には、収束条件を満足したと判断し、ステップＳ２０５で、上記ステップＳ２０３で更新された各時刻ｋのフレーズ指令ｕp［ｋ］及びアクセント指令ｕa［ｋ］と、上記ステップＳ２０１で更新された事後確率Ｐ(sk=t|y,o′,θ′)とに基づいて、上記式（７３）、式（７４）、又は式（７５）、式（７６）に従って、各時刻ｋのフレーズ指令の振幅A(p)[k]、及び各位置ｎのアクセント指令の振幅Aa(a)を更新することにより、各時刻における状態に応じたフレーズ指令の振幅及び各アクセント指令の振幅を表すパラメータθを更新する。
【０２７３】
ステップＳ２０６において、予め定められた収束判定条件を満たしたか否かを判定し、収束判定条件を満たしていない場合には、上記ステップＳ２０１へ戻る。一方、収束判定条件を満たした場合には、ステップＳ２０７において、
【０２７４】
上記ステップＳ２０３で最終的に更新された指令系列ｏに基づいて、Viterbi アルゴリズムにより、状態系列ｓを推定する。また、推定された状態系列ｓから、状態遷移確率φを求める。
【０２７５】
そして、ステップＳ２０８において、上記ステップＳ１０３で学習されたパラメータγ、上記ステップＳ２０５で最終的に得られたパラメータθと、上記ステップＳ１０６で得られた状態遷移確率φとを、パラメータ記憶部４０に格納する。
【０２７６】
次に、予測対象のソース音声の時系列データが、基本周波数パターン予測装置１００に入力されると、基本周波数パターン予測装置１００において、図１２に示す基本周波数パターン予測処理ルーチンが実行される。
【０２７７】
まず、ステップＳ１２１において、入力された予測対象のソース音声の時系列データを読み込み、各時刻ｋのスペクトル特徴量ベクトルｃ［ｋ］を抽出する。ステップＳ１２２において、パラメータ記憶部４０に記憶されたパラメータγと、上記ステップＳ１２１で抽出された各時刻のスペクトル特徴量ベクトルc[k]とに基づいて、各時刻の基本周波数ｙ［ｋ］を初期設定する。
【０２７８】
そして、ステップＳ１２３では、パラメータ記憶部４０に記憶されたパラメータγと、上記ステップＳ１２１で抽出された各時刻のスペクトル特徴量ベクトルc[k]とに基づいて、各時刻ｋの正規分布のインデックス＾ｍkを推定する。
【０２７９】
ステップＳ２２１では、上記ステップＳ１２２で初期設定された、または後述するステップＳ１２５で前回更新された各時刻ｋの基本周波数ｙ［ｋ］を固定して、上記ステップＳ２０１〜ステップＳ２０６と同様に、事後確率Ｐ（ｏ，θ｜ｙ）を局所最大化する指令系列ｏを推定する。
【０２８０】
そして、ステップＳ２２２において、上記ステップＳ２２１で推定された指令系列ｏと、上記ステップＳ１２３で推定された各時刻の正規分布のインデックス＾ｍkとに基づいて、上記式（１０１）に従って、各時刻ｋの基本周波数ｙ［ｋ］を更新する。
【０２８１】
ステップＳ１２６では、予め定められた収束判定条件を満たしたか否かを判定し、収束判定条件を満たしていない場合には、ステップＳ２２１へ戻る。一方、収束判定条件を満たした場合には、ステップＳ１２７において、上記ステップＳ２２２で最終的に得られた各時刻ｋの基本周波数ｙ［ｋ］を、予測対象のソース音声に対応するターゲット音声の各時刻の基本周波数の予測結果として、出力部９０により出力し、基本周波数パターン予測処理ルーチンを終了する。
【０２８２】
以上説明したように、第２の実施の形態に係る基本周波数パターン予測装置によれば、ソース音声の各時刻のスペクトル特徴量ベクトルと、ターゲット音声の各時刻の基本周波数との間の関係をモデル化したＧＭＭである第１確率分布Ｐ(q[k]，c[k]｜γ)のパラメータγを学習し、基本周波数パターン生成過程をモデル化した第２確率分布Ｐ（ｙ，ｏ｜θ，φ）のパラメータθ、φを学習し、予測対象のソース音声の時系列データから抽出される各時刻のスペクトル特徴量ベクトルから、第１確率分布Ｐ(q[k]，c[k]｜γ)と第２確率分布Ｐ（ｙ，ｏ｜θ，φ）との積を用いて表される規準を大きくするように、予測対象のソース音声に対応するターゲット音声の各時刻の基本周波数を予測することにより、F0 パターンの物理的な生成過程の制約を考慮しながらスペクトル特徴量系列に対応する最適なF0パターンを推定することができる。
【０２８３】
＜実験＞
【０２８４】
図１３に示すF0パターンの音声データに対し、上述した従来手法である関連技術１のスペクトル特徴量系列からのF0 パターン予測方法と、本発明の第１の実施の形態に係る手法とによりスペクトル特徴量系列からF0 パターンの予測を行う実験を行った。図１４に、両手法により予測されたF0パターンを示す。図１４では、実線が、従来手法による音声特徴量系列からのF0 パターンの予測結果の例を示し、点線が、第１の実施の形態に係る手法による音声特徴量系列からのF0 パターンの予測結果の例を示す。
【０２８５】
図１３のF0パターンとの近さがF0パターンの良さの指標になる。そこで、それぞれの手法で得られたF0パターンと、図１３のF0パターンとのコサイン距離（1 に近いほど近いこと意味する）を測ったところ、従来手法が0.55、第１の実施の形態に係る手法が0.59 であった。このことから、本発明の第１の実施の形態の手法の、従来手法に対する優位性が示された。
【０２８６】
なお、本発明は、上述した実施形態に限定されるものではなく、この発明の要旨を逸脱しない範囲内で様々な変形や応用が可能である。
【０２８７】
例えば、上述の基本周波数パターン予測装置は、内部にコンピュータシステムを有しているが、「コンピュータシステム」は、ＷＷＷシステムを利用している場合であれば、ホームページ提供環境（あるいは表示環境）も含むものとする。
【０２８８】
また、本願明細書中において、プログラムが予めインストールされている実施形態として説明したが、当該プログラムを、コンピュータ読み取り可能な記録媒体に格納して提供することも可能である。
【符号の説明】
【０２８９】
１０   入力部
２０   演算部
３０   学習部
３２   特徴量抽出部
３４   基本周波数系列抽出部
３６   第１モデルパラメータ学習部
３８   第２モデルパラメータ学習部
４０   パラメータ記憶部
５０   変換処理部
５２   特徴量抽出部
５４   基本周波数系列予測部
５６   正規分布系列予測部
５８   状態系列推定部
６０   基本周波数系列更新部
６２   収束判定部
９０   出力部
１００ 基本周波数パターン予測装置
２５６ 正規分布系列予測部
２５８ 指令系列推定部
２６０ 基本周波数系列更新部
【図１】
0006468519.tif 000085
【図２】
0006468519.tif 000086
【図３】
0006468519.tif 000087
【図４】
0006468519.tif 000088
【図５】
0006468519.tif 000089
【図６】
0006468519.tif 000090
【図７】
0006468519.tif 000091
【図８】
0006468519.tif 000092
【図９】
0006468519.tif 000093
【図１０】
0006468519.tif 000094
【図１１】
0006468519.tif 000095
【図１２】
0006468519.tif 000096
【図１３】
0006468519.tif 000097
【図１４】
0006468519.tif 000098