2-11-6
ノート指令と表現指令によって駆動される歌声 F0 生成過程の統計モデル∗◎大石康智 (NTT), 亀岡弘和 (NTT/東大), 持橋大地 (統数研), 柏野邦夫 (NTT)1

はじめに歌声の音高 (F0 ) 軌跡には，発声器官の物理的制約に起因する動的変動成分（オーバーシュートや微細変動成分）や歌唱者の意図的表現による動的変動成分（ビブラートやポルタメント）が含まれ，知覚的には前者は人間らしさ・自然性に関係し，後者は巧拙感に関係することがわかってきている [1]。ただし，これらの動的変動成分を F0 軌跡から特徴抽出し，歌唱者ごとにその特性を精緻に学習することはまだ十分に検討されていない。
本研究の目的は，このような動的変動成分を F0 軌跡から楽譜情報と分離して特徴抽出し，歌唱者ごとにどのようなパターンをもちうるのか，各パターンが文脈（楽譜の音符列）にどう依存するかを計算機に学習させることである。ここでは，うろ覚えの状態で歌った歌声ではなく，楽曲のメロディまたはその楽譜を既知として，歌唱者なりに表情付けして歌った歌声を対象とする。本研究は，歌唱者の歌い方や個性，癖を学習することを目指しており，歌唱力評価や歌唱者識別，現在盛んに研究される歌声合成や歌声変換[2, 3, 4] への応用が期待できる。
本研究では，物理モデルに基づいて F0 軌跡の生成過程を記述し，そこから動的変動成分の特徴抽出に取り組む。これまで，話声の F0 パターンを表現する藤崎モデル [5] を参考に，2 次系を利用した歌声 F0 制御モデルが提案された [1, 6]。2 次系の伝達関数はG(s) =Ω2s2 + 2ζΩs + Ω2(1)であり，この減衰率 ζ によって，様々な振動現象を表現できる。文献 [1] では，楽譜の音符列を表す階段状信号に式 (1) のインパルス応答を部分的に畳み込んで得られる F0 軌跡を利用して，表情豊かな歌声合成音を実現した。制御パラメータ（減衰率 ζ と固有周波数 Ω）は手作業あるいは規則に基づいて決定された。
一方，我々は観測される F0 軌跡から制御パラメータを推定する逆問題の解法を検討してきた [7]。ただし，ビブラートや微細変動成分がすべて白色雑音としてモデル化されたため，意図的表現による動的変動成分を微細変動成分と分離して特徴付けられなかった。
本稿では，歌唱の意図的表現を特徴抽出するために，ノート指令信号と表現指令信号によって駆動される歌声 F0 軌跡の生成過程を提案する（Fig. 1）。ここ
で，ノート指令信号は楽譜に記載される音符の並びを∗Fig. 1

提案する歌声の F0 生成過程の概略表現する。一方，表現指令信号は歌唱者の音楽的な表現意図を矩形状の細かい指令として表現する。ノート成分と表現成分は各指令信号によって駆動される制御機講の出力であり，これらの制御機講は 2 次系で表現される。ノート制御機講はオーバーシュートなどのノートの立ち上がり方を制御する。表現制御機講は矩形状の細かい指令信号を制御してビブラートやポルタメントを生成する。微細変動成分は 10 Hz 以上の不規則な振動成分を想定する [1]。最終的に，対数スケールの F0 軌跡 y(t)（ここで，t は時間）は，これら 3 つの成分の重ね合わせと想定する。この F0 生成過程を想定した理由は 2 点ある。
理由 1 藤崎モデルでは，甲状軟骨の二つの独立な運動 (平行移動と回転) に伴う声帯の長さの変化の合計が F0 の時間的変化をもたらすと解釈され，これらの運動に関係する成分をそれぞれ，フレーズ成分とアクセント成分とした [5]。歌声にも同様の機講が存在すると仮定し，大域的な変化を表すノート成分とアクセントのような歌唱者が意図的に制御できる表現成分から歌声 F0 軌跡が構成されると考えた。
理由 2 従来，ビブラートは正弦波で表現され，そのパラメータは Vibrato rate と Vibrato extent とされた [1]。これらのパラメータが指数的に変化するモデルも提案された [3]。しかし実際は，これらのパラメー
タは時々刻々と変化するものの，指数的な変化であるとは言えない。この複雑な振動現象をモデル化して歌い方を認識するために，意図を表す矩形状の細かい指令信号が 2 次系によって制御されるビブラートの生成過程を想定した。
以降の節では，上記の歌声 F0 軌跡の生成過程を離散時間表現して確率モデル化し，統計的手法に基づいてモデルパラメータの推定アルゴリズムを導出する。そして，提案モデルの有効性を実測の F0 軌跡を用いて評価する。文献 [8] では，同じ枠組みで藤崎モデルを確率過程に基づいて統計モデル化する試みも行っているのでこちらも参照されたい。

A Statistical Model of Singing Voice F0 Contours Driven by Musical Note and Expression Command Functions by OHISHI, Yasunori (NTT), KAMEOKA, Hirokazu (NTT / The University of Tokyo), MOCHIHASHI, Daichi (The Institute of Statistical Mathematics), KASHINO Kunio (NTT)

日本音響学会講演論文集- 345 -2012年3月

2

歌声 F0 生成過程の確率モデル

歌声 F0 軌跡の生成過程を離散時間表現し，確率過程に基づいて統計モデル化する。連続時間領域で表現される 2 次系制御機講の伝達関数の離散時間表現を得るために，後退差分変換を利用する [7]。t0 をサンプリング周期とすると，この変換によりノート制御機構の逆システムの伝達関数 Hn−1 (s) は z 領域で，Hn−1 (z) = a2 z −2 + a1 z −1 + a0(2)と書くことができる。ただし，a2 = ϕ2 , a1 = −2ϕ(ψ+ϕ), a0 = 1 + 2ϕψ + ϕ2 ，および，ϕ = 1/(Ωt0 ), ψ = ζと表現される。ここで，k を離散時刻インデックスとし，ノート指令信号およびノート成分の離散時間表
現をそれぞれ un [k]，yn [k] とすると，yn [k] は，ノート制御パラメータ ϕ, ψ によって特性が決まる拘束つき全極モデルからの出力

un [k] = a0 yn [k] + a1 yn [k − 1] + a2 yn [k − 2](3)と見なすことができる。同様に，表現指令信号 ue [k]と表現成分 ye [k] の関係も ue [k] = b0 ye [k] + b1 ye [k −1] + b2 ye [k − 2] と書ける。ただし，表現制御機講は臨界制動系 (ζ = 1 の場合) を想定し，b2 = ξ 2 , b1 =−2ξ(1 + ξ), b0 = 1 + 2ξ + ξ 2 ，表現制御パラメータを ξ = 1/(Ωt0 ) とする。微細変動成分 yf (t) の離散時間表現を yf [k] として，歌声 F0 軌跡の離散時間表現
は，y[k] = yn [k] + ye [k] + yf [k] と想定する。
次に，F0 生成過程を確率モデル化する。ノート指令信号と表現指令信号はそれぞれ，楽譜に記載されるメロディの音符の並びと歌唱者の音楽的な表現意
図を表す (Fig. 1）。これらの指令信号を表現するために，HMM を利用して，un [k] と ue [k] を確率モデル化する（Fig. 2）。o[k] := (un [k], ue [k])T を"#"#µn [k]
σn2 0, Υ :=o[k] ∼ N (ν[k], Υ), ν[k] :=µe [k]
0 σe2のように正規分布する確率変数と見なし，平均 ν[k]が状態遷移に伴って変化するモデルを考える。
具体的には，この HMM は I × J 個の状態集合 S :=(i){Si,j }I,Ji=1,j=1 からなる。状態 Si,j では，µn [k] は An +(i)di の値をとる。ここで，An は楽譜に記載されるメロ
ディの i 番目の音符の音高を表し（この値は楽譜から与えられるものとする），di はその絶対音高からのずれ（音高シフトパラメータと呼ぶ），I は音符の総数を表す。Fig. 2 に示す状態遷移により，µn [k] は I 個の音符区間からなる階段状信号を表現する。一方，µe [k] は(i,j)
Be の値をとり，これは歌唱者の表現意図を表すための，矩形の指令の大きさを表す。ここで，i 番目の音符では，Si,1 を通らずして，状態 Si,j から別の状態 Si,j 0(j 6= j 0 , 2 ≤ j ≤ J, 2 ≤ j 0 ≤ J) へ直接に遷移できない制約を設けることによって，µe [k] は Fig. 2 に示すよ日本音響学会講演論文集

Fig. 2HMM による指令信号の統計モデル化

うな矩形信号を表現する。J は表現指令を構成するための状態数を表す（Fig. 2 は J = 3 である）。指令信号を生成する HMM（指令信号モデル）の構成をまとめI,Jると，出力系列: {o[k]}K
k=1 ，状態集合: {Si,j }i=1,j=1 ，K
状態系列: {sk }k=1 ，出力確率分布: P (o[k]|sk = n) =(i)(i,j)
N (cn [k], Υ), cn [k] = (An +di , Be )T , (n = Si,j )，状態遷移確率 φn0 ,n := log P (sk = n|sk−1 = n0 ) となる。簡単のため，状態遷移確率 φn0 ,n は定数とすると，指令信号モデルにおいて推定すべきパラメータ
は，HMM の経路を表す状態系列 {sk }K
k=1 ，音高シフI
トパラメータ {di }i=1 ，表現指令の大きさパラメータ(i,j)
22{Be }I,Ji=1,j=1 ，指令信号の分散パラメータ σn , σeであり，これらをまとめて θu と記述する。また，平K
均値系列 {µn [k]}K
k=1 および {µe [k]}k=1 は，状態系
T列 {sk }K
k=1 が与えられたもとで，(µn [k], µe [k]) ←csk [k] で与えられる。
y = (y[1], . . . , y[K])T の確率密度関数を導く。un:=(un [1], . . . , un [K])T , ue:=(ue [1], . . . , ue [K])T , µn :=(µn [1], . . . , µn [K])T , µe := (µe [1], . . . , µe [K])T とするとun |θu ∼ N (µn , Σn ),Σn = σn2 I K(4)ue |θu ∼ N (µe , Σe ),Σe = σe2 I K(5)と書ける。ここで，I K は K × K の単位行列を表す。ノート成分 y n := (yn [1], . . . , yn [K])T とノー
ト 指 令 信 号 un の 関 係 ，お よ び 表 現 成 分 y e :=(ye [1], . . . , ye [K])T と表現指令信号 ue の関係は，
 (1)   (1)a0OA(1) (1)
 .  a1 a0
 ..   .
.
.
.
.
.
.
.
.
A :=  A(i)  = (i) (i) (i)a2 a1 a0 ..  .
.
.
 .  .. .. ..
(I) (I) (I)A(I)
Oa2 a1 a0
b0O
b1 b0bbb
B :=  2 1 0
 ... ... ... 
Ob2 b1 b0- 346 -

2012年3月

と置くと，それぞれ，un = Ay n , ue = By e と表現できる。ここでは，先行研究で得られた知見 [1] に基づいて，ノート制御パラメータがノートごとに異なるも(i)(i)(i)(i)(i)のと想定し，a0 , a1 , a2 は，a2 = (ϕ(i) )2 , a1 =(i)−2ϕ(i) (ψ (i) + ϕ(i) ), a0 = 1 + 2ϕ(i) ψ (i) + (ϕ(i) )2 と(i)(i) I表現され，{ϕ , ψ }i=1 を制御パラメータとする。
そこで，ノート成分と表現成分の確率密度関数は，

y n ∼ N (A−1 µn , A−1 Σn (A−1 )T )(6)

y e ∼ N (B −1 µe , B −1 Σe (B −1 )T )(7)が導出される。微細変動成分 y f := (yf [1], . . . , yf [K])Tはガウス性白色雑音を想定する。

y f |σf2 ∼ N (0, Σf ),Σf = σf2 I K(8)仮定より，y n , y e , y f は独立なので，パラメータΘ := {θu , {ϕ(i) , ψ (i) }Ii=1 , ξ, σf2 } が与えられた下での歌声 F0 軌跡 y = y n + y e + y f の確率密度関数（尤度関数）は式 (6)，(7) と式 (8) より，¡y|Θ ∼ N A−1 µn + B −1 µe ,(9)¢−1−1 T
−1−1 TA Σn (A ) + B Σe (B ) + Σfとなる。Θ の事前確率は，各パラメータは互いに独立で，パラメータ σn2 , σe2 , σf2 は一様に分布すると仮定し，P (Θ) ∝ P (ϕ(1) ) . . . P (ϕ(I) )P (ψ (1) ) . . .
(I,J)
(1,1)
)
P (ψ (I) )P (ξ)P (d1 ) . . . P (dI )P (Be ) . . . P (BeQK
P (s1 ) k=2 P (sk |sk−1 ) と記述される。

3

パラメータ推定アルゴリズム

Θ の事後確率 P (Θ|y) を最大化する問題を解析的T
T Tに解くことはできないが，x := (y Tn , ye , yf ) を完全データと見なすことで，EM アルゴリズムによる不完全データ問題に帰着できる。不完全データと完全データの関係は y = Hx となり，ここで H :=[I K I K I K ] とする。このとき，Q 関数は，
1£log |Λ−1 | − tr(Λ−1 E[xxT |y; Θ0 ])2¤+ 2mT Λ−1 E[x|y; Θ0 ] − mT Λ−1 m + log P (Θ)
 −1  T −1A µnA Σn AOOm := B −1 µe  , Λ−1:=OB T Σ−1e B O 0OOΣ−1f
Q(Θ, Θ0 ) =cとなる。ここで，tr(·) は行列のトレースを表す。E ス
テップでは，直前のステップで更新されたパラメー
タを Θ0 に代入し，条件付きガウス分布の性質から
E[x|y; Θ0 ] と E[xxT |y; Θ0 ] を計算する [7]。y n , y e , y fに対応するように，E[x|y; Θ0 ] および E[xxT |y; Θ0 ] を以下のように区分表現する。
 
x̄nRn ∗ ∗ E[x|y; Θ0 ] =  x̄e  , E[xxT |y; Θ0 ] =  ∗ Re ∗ x̄f∗ ∗ Rf

日本音響学会講演論文集各パラメータの M ステップ更新式を以下に示す。
1) 状態系列: Q 関数で s := {sk }K
k=1 に関する項はK

I1 (s) := −1X(o[k] − csk [k])T Υ−1 (o[k] − csk [k])2k=1+ log P (s1 ) +K
Xlog P (sk |sk−1 ) (10)k=2である。ただし，o[k] := ([Ax̄n ]k , [B x̄e ]k )T であり，[·]k はベクトルの k 番目の要素を表す。これを最大化する {sk }K
k=1 は Viterbi 探索により導出できる。
2) ノート制御パラメータ: ϕ(i) と ψ (i) に関する事
2
前分布を ϕ(i) ∼ N (µϕ , σϕ)，ψ (i) ∼ N (µψ , σψ2 ) とす(i)る。Q 関数の中で ϕ と ψ (i) に関する項は，

I2 (ϕ(i) , ψ (i) ) = |TSi,· | log(1 + 2ϕ(i) ψ (i) + (ϕ(i) )2 )1
1− 2 tr((A(i) )T A(i) Rn ) + 2 ([µn ]TSi,· )T A(i) x̄n2σnσn1
1− 2 (ϕ(i) − µϕ )2 − 2 (ψ (i) − µψ )22σϕ2σψTSi,· = {k|sk ∈ {Si,1 , . . . , Si,J }} (11)となる。ここで，|T | は集合 T の要素数を表す。また，[µ]T は，集合 T の要素を添え字として，その添え字に相当する µ の要素を取り出した部分ベクトルを表す。今，定数行列 U 0 , U 1 , U 2 を用いて A(i) は，A(i) = (ϕ(i) )2 [U 2 ]TSi,· +ϕ(i) ψ (i) [U 1 ]TSi,· +[U 0 ]TSi,·と表現される。ここで，[U ]T は集合 T の要素を添え字として，行列 U からその添え字に相当する行ベクトルを取り出して構成される部分行列を意味する。
ニュートン・ラフソン法を利用して，I2 (ϕ(i) , ψ (i) ) を最大化する ϕ(i) と ψ (i) が数値的に導出される。
3) 表現制御パラメータ: ノート制御パラメータと同様の手順で更新式を導出できる。ただし，推定すべきパラメータ数が多いため，今回は ξ = 3 に固定する。
(i,j)4) その他のパラメータ: di と Be の事前分布をそ(i,j)
22∼ N (µB (i,j) , σB) とす
れぞれ di ∼ N (0, σd ) と Beる。残されたパラメータの更新式は
X
1
di =([Ax̄n ]k − A(i)(12)n )
22|TSi,· | + σn /σdk∈TSi,·

Be(i,j) =³ X [B x̄ ]1µB (i,j) ´e k+
22|TSi,j |/σe2 + 1/σB
σe2σBk∈TSi,jTSi,j = {k|sk = Si,j }(13)´
Tσn2 = tr(AT ARn ) − 2µTn Ax̄n + µn µn /K (14)´³T
σe2 = tr(B T BRe ) − 2µT
Bx̄
+µµeee e /K (15)³σf2 = tr(Rf )/K(16)と導出される。無声音のため F0 が推定されない区間の取り扱いや各パラメータの実際の更新方法については，文献 [9] を参照されたい。

- 347 -

2012年3月

4

評価実験

実測の F0 軌跡を用いて，提案手法の性能を評価する。「RWC 研究用音楽データベース：ポピュラー音楽」(RWC-MDB-P-2001)[10] における，歌手名：緒方智美，曲番号：No.07，曲名：PROLOGUE のメロ
ディの F0 を手作業でラベル付した結果 [11] を実測 F0軌跡として利用した。本来ならば音響信号から F0 を推定すべきであるが，今回は提案手法の性能の上限を調べるためにこのようなデータを利用した。F0 は10ms ごとにラベル付されているので，アップサンプリングによって，5ms ごとの実測 F0 軌跡を得た。同楽曲の楽譜情報はその MIDI データを利用した。
パラメータ推定アルゴリズムの設定値を以下に示(i)す。I と {An }Ii=1 は MIDI データから与えられる。
表現指令信号を構成するための状態数 J は 5 とした。
HMM の状態遷移確率は，φSi,1 ,Si,1 = log(0.9999 ×(J − 1)/J), φSi,1 ,Si,j = log(0.9999/J), φSi,1 ,Si+1,1 =log(0.0001), φSi,j ,Si,j = log(0.9999), φSi,j ,Si,1 =log(0.0001), (1 ≤ i ≤ I, 2 ≤ j ≤ J) とした。
パラメータの事前分布における固定パラメータは，µϕ = 6, σϕ2 = 0.1, µψ = 0.6, σψ2 = 0.02, σd2 =22500, σB= 100, µB (i,1) = 0, µB (i,2) = 30, µB (i,3) =−30, µB (i,4) = 60, µB (i,5) = −60, (1 ≤ i ≤ I) とした。これらは文献 [1] と予備実験に基づいて決定した。
実測 F0 軌跡からノート指令信号と表現指令信号，制御機構の出力であるノート成分と表現成分の推定結果を Fig. 3 に示す。これらを足し合わせたものを合成 F0 軌跡とし，定性的には実測 F0 軌跡に近い軌跡が得られた。また，微細変動成分パラメータである
σf2 の値は 187.7 であるため，振幅が 13 cent 程度の変動成分が推定され，先行研究の知見 [1] と整合する。
無声音の区間を欠損データとして扱うため，合成 F0軌跡ではその部分が周囲の F0 から補間される。
Fig. 4 は具体的に，ポルタメントとビブラートが生成される様子を示す。楽譜（MIDI）では一つの音を伸ばし続けるよう記述されているものの，歌手の表現意図によって，半音だけ滑らかに変化させるポルタメントが観測される。このとき，ノート指令信号は同一音符区間として推定され，表現指令信号は音高を段階的に変化させる指令が推定される。一方，実測 F0軌跡には必ずしも正弦波とは言えないビブラートが観測される。提案法では，部分的に矩形変化する表現指令が推定されてビブラートを生成する。ポルタメ
ントとビブラートの生成のどちらにおいても合成 F0軌跡は実測 F0 軌跡に近い軌跡が得られた。

5

おわりに歌声 F0 軌跡に含まれる動的変動成分を楽譜情報と分離して抽出することを目的として，2 つの指令信号によって駆動される F0 軌跡の生成過程を記述し，モ

日本音響学会講演論文集

Fig. 3 実測 F0 軌跡に対するパラメータ推定結果

Fig. 4

ポルタメントとビブラートの生成デルパラメータの推定アルゴリズムを提案した。評価実験では，実測の F0 軌跡から推定される歌唱者の表現意図を考察した。提案モデルの有効性を主観的に評価する聴取実験も行なっているので，文献 [9] も参照されたい。今後の課題は，抽出された動的変動成分のパターンを歌唱者ごとに学習することである。

参考文献[1] Saitou et al., Proc. EUROSPEECH’09, pp.
832–835, 2009.
[2] Mase et al., Proc. INTERSPEECH’10, pp.
845–848, 2010.
[3] 右田ほか, 情処論文誌, Vol.52, No.5, pp. 1910–1922, 2011.
[4] 中野ほか, 情処論文誌, Vol.52, No.12, pp. 3853–3867, 2011.
[5] Fujisaki, In Vocal Physiology: Voice Production, Mechanisms and Functions, (O. Fujimura,ed.) Raven Press, pp. 347–355, 1988.
[6] 柏野ほか, 音講論集, 2-9-1，pp. 625–626 1998.
[7] 大石ほか, 音講論集, 1-8-19 pp. 279–282, 2011.
[8] Kameoka et al., Proc. SAPA’10, pp. 43–48,2010.
[9] 大石ほか, 情処研報, 2012–MUS94–12, 2012.
[10] 後藤ほか, 情処論文誌, Vol.45, No.3, pp. 728–738, 2004.
[11] Goto, Proc. ISMIR’06, 2006.

- 348 -2012年3月