Vol.2010-MUS-86 No.7
2010/7/28

情報処理学会研究報告
IPSJ SIG Technical Report

1. は じ め に

混合音中の歌声スペクトル包絡推定に基づく歌声の声質変換手法

本稿では，混合音中の歌声の声質の変換手法について述べる．つまり，入力として伴奏を含む歌声の音響信号と変換先歌手の歌声の音響信号を取り，歌声の声質が変換された音響信号を出力する手法である．ここで声質とは，歌声のスペクトルの静的な形状のことを指し，

藤 原 弘

将†1

後

藤

真

孝†1

基本周波数（F0）の動きなど，動的な成分は含まないものとする．変換先歌手の歌声の音
響信号は、変換元のものと同じ楽曲で有る必要はなく、複数の楽曲でも良い。一方で、単独歌唱の音響信号である必要があり、また変換元の音響信号に含まれる歌声の母音が含まれて本稿では，混合音中の歌声の声質を，別の歌手の声質に変換する手法について述べる．従来の声質変換は単独歌唱のみを対象としていた。本研究では，藤原ら1) によって提案された W-PST 法を応用する．W-PST 法によりスペクトル中の歌声が優勢な
周波数成分を同定できる。本稿では、まず混合音のスペクトルに対して、W-PST 法で推定された歌声が優勢な周波数成分を操作することで、歌手の声質を他の歌手のものに変換する手法について述べる。次に、W-PST 法の実行に必要なスペクトル包絡を、伴奏が混在した歌声から推定する手法について述べる。本手法を実装し，実際の楽曲に対して適用したところ，歌声の声質が変換できることを確認した．

いる必要がある。
近年，能動的音楽鑑賞インタフェース2) と呼ばれる，音楽を自分好みに操作しながら，より能動的に音楽を鑑賞するための技術とインタフェースが提案されている．例えば，吉井らによる Drumix3) は，ドラムの音量調整とドラムのパターンを置き換えが可能で，糸山らによる Instrument Equalizer4) では，各パートごとの音量を自由に操作しながら可能であった．本研究の技術はそのような能動的音楽鑑賞を歌声に対して実現する技術として位置づけることができ，楽曲の歌手の声質を自分好みの歌手の声に置き換えて鑑賞することが可能に

Singing voice conversion methodby using spectral envelope of singing voiceestimated from polyphonic music

なる．
また，歌声合成技術の発達や Web 上の動画共有サイトの発達により，一般ユーザーが音楽を作成するようになってきており，それを支援するための製品や技術が登場している．例
として，Bonada らの素片連結型歌声合成技術5) に基づく YAMAHA 社の歌声合成ソフト

Hiromasa Fujihara†1 and Masataka Goto†1

ウェア VOCALOID6) や、酒向らの HMM 歌声合成技術7) に基づく歌声合成ウェブサービ
ス Sinsy8) などがある。このような技術により，ユーザーは多くの人とコラボレーションしThis paper describes a singing voice conversion method that can deal withsinging voices in polyphonic music. Conventional voice conversion methods
only deal with monophonic singing voices. In this paper, we utilize the W-PST
method proposed by Fujihara et al.1) , which can identify the frequency components of a singing voice in a polyphonic spectrum. We first describe our methodof converting the vocal timbres of singing voices to those of other singers bymaniplulating only frequency components of singing voices identified by the
W-PST method. Since the W-PST method requires spectral envelopes of thesinging voices, we then describe a method of estimating them from polyphonicmusic. We applied our method to actual musical audio signals and confirmedthat it was able to convert the vocal timbre of the singing voices in polyphonicmusic.

ながら、より手軽に楽曲を制作し、作品を発表することができるようになっている．本研究の技術は，ユーザーが既存の伴奏を含む楽曲の声質を別の歌手や自分の声に変えることができるため，新たな音楽制作ツールとして使用できる可能性も秘めている．
他の音が背景等として含まれないクリーンな話し声を対象とした声質変換は数多く研究がなされており9)–11) ，これらの技術の一部は単独歌唱の歌声にも適用が可能である．また、河原らの開発した歌声分析合成システム STRAIGHT に基づく歌声のモーフィング12) では，

2 種類の単独歌唱の歌声をリアルタイムにモーフィングし，ある歌手の声質で別の歌手の歌†1 産業技術総合研究所
National Institute of Advanced Industrial Science and Technology (AIST)1

c 2010 Information Processing Society of JapanVol.2010-MUS-86 No.7
2010/7/28

4000

歌声スペクトルテンプレート Yv, fLogPower

合成を表す

LogPower

⊕

500
1000
2000
Frequency (Hz)

図1

(a)

4000

ノイズスペクトルテンプレート Y

n, f

500
1000
2000
Frequency (Hz)4000(c) 確率的スペクトルテンプ
レート Y

500
1000
2000
Frequency (Hz)

歌声包絡テンプレート Y '(d)LogPower

4000
v, f

F0の値を調整して、様々なピッチの音を表現

LogPower

+

図2

gn
500
1000
2000
Frequency (Hz)(b)

f0
LogPower

gv
500
1000
2000
Frequency (Hz)(a)

重みパラメータを調整し様々なSIRのス
ペクトルを表現
歌声スペクトルテンプレートとノイズス
ペクトルテンプレートを合成して確率的スペクトルテンプレートができる
確率的スペクトルテンプレートから観測
スペクトルが生成される
LogPower

LogPower

情報処理学会研究報告
IPSJ SIG Technical Report

(b)

500
1000
2000
Frequency (Hz)4000

駆動音源関数 H ( f ; f )
0

(c)

500
1000
2000
Frequency (Hz)4000

歌声スペクトルテンプレート Yv, f

歌声スペクトルテンプレートの例1) ．歌声包絡テンプレートと駆動音源関数から生成される．

2. W-PST 法に基づく混合音中の歌声の声質変換

4000

観測スペクトル y( f )本研究では，混合音中の歌声の声質変換を実現するために，藤原ら1) によって提案された

f

W-PST 法を応用する．W-PST 法とは，連続ウェーブレット変換（CWT）によって得られた観測スペクトルを歌声と歌声以外の音（ノイズ）が混ざった状態としてモデル化し，音観測スペクトルの生成過程1) ．図の濃淡は確率密度を表現する．重みパラメータ gv と gn を調整することで，様々な SIR のスペクトルを表現できる．

素と F0 を推定する手法である．本研究は，W-PST 法は歌声とノイズの SIR（Signal-to-

Interference Ratio）を推定するため，混合音のスペクトル中で歌声の周波数成分が優勢ない方の歌を作成することなどができる．また、この STRAIGHT のモーフィングを声質変換に応用した研究例もある

13)周波数帯域を同定できることを利用し，混合音中の歌声の声質変換に応用する．本手法は，．しかし，これらの単独歌唱を対象とした技術は混合音には変換元の音響信号と変換先の歌手の音響信号を入力とし，変換元の音響信号の歌声の声質を適用できず，伴奏を含む混合音中の歌声の声質変換は今まで扱われてこなかった．

変換先の歌手のものに変換した音響信号を出力する．本稿では，変換元および変換先の音響

混合音中の歌唱の声質変換を行う際の本質的難しさは，歌声を処理する際に伴奏音の影響

信号について，音素と F0 のラベル（各時刻における音素名と F0 の値）が付与されているを排除する必要があるだけではなく，歌声への処理が伴奏音に与える影響を排除する必要がことを仮定し，変換先の音響信号は単独歌唱のものであると仮定する．ただし、音素と F0ある点である．なぜなら，歌声以外の音の音質を劣化させずに，歌声のみの音質を変化させのラベルは文献 1) の手法で推定することが可能であり、自動推定したラベルを使用して声1)によって提案された W- PST 法を応

質変換を行うことに今後取り組む予定である。また，本章の以下の処理は母音区間に対して用して，歌声の周波数成分のみを操作することを可能にした．W-PST 法は，混合音中の歌のみ行われる。実際は、子音にも個人性が存在するため、子音区間に対して処理を行うこと

声の F0 と音素を推定する手法で，伴奏音と歌声が混ざった状態としてモデル化し，歌声のは今後の課題である。

る必要があるからである．そこで本研究は，藤原ら

2.1 概周波数成分が優勢な帯域を同定することが可能である．しかし，W-PST 法では歌声のスペ

W-PST 法

クトル包絡推定は，単独歌唱のデータを用いていた．一方，声質変換の目的では，変換した

要1)では，図 1 (c) と (d) で示されるように，歌声を含む混合音のスペクトルがい音響信号は伴奏が混在した混合音として与えたいことが多いので，そのままではスペクト

確率的スペクトルテンプレートと呼ばれる確率分布の集合から生成されると仮定する．さル包絡推定ができない．そこで，本研究では，混合音の音響信号中の歌声のスペクトル包絡らに，パワースペクトルの加法性を仮定し，確率的スペクトルテンプレートを，歌声を表推定手法を新たに開発することで，この問題を解決した．

現する歌声スペクトルテンプレート（図 1 (a)）と歌声以外の音を表現するノイズスペクト
ルテンプレート（図 1 (b)）の加算で表現する．つまり，観測スペクトルを生成する音源を，

2

c 2010 Information Processing Society of JapanVol.2010-MUS-86 No.7
2010/7/28

情報処理学会研究報告
IPSJ SIG Technical Report

換先音響信号を取る。変換元音響信号は、市販 CD 等の実世界の音楽音響信号で、通常歌変換先音響信号（単独歌唱）

変換元音響信号（混合音）

声と共に伴奏音が含まれている。一方、変換先音響信号は声質を変える目標となる歌手の音響信号で、本研究では単独歌唱であることを仮定する。前処理としてこれらの音響信号から

CWT によりスペクトログラムを計算しておく（図 3 (a), (b)）。本研究では、まず W-PST(b) 連続ウェーブレット変換(a) 連続ウェーブレット変換

法により変換元の観測スペクトルを歌声スペクトルテンプレートとノイズスペクトルテン

変換先観測スペクトル(c) 変換元歌手のテンプレートを推定

プレートの合成で表現することを目指す。これにより、変換元スペクトルの周波数成分ごとに、歌声が優勢なのか、伴奏成分が優勢なのかを同定することができる。そのためには、変換元観測スペクトル(d) SIRの推定
変換元歌声テンプレート(e) 歌声テンプレートを変換先
歌手のものと置き換える

まず変換元の観測スペクトルを表すのに最適な歌声包絡テンプレートとノイズスペクトル(f) 変換先歌手のテンプレートを推定

テンプレートを変換元の観測スペクトル自体から推定する必要がある（図 3 (c)）。なぜな
ら、変換元スペクトル中の歌手の声質は事前に単独歌唱の学習データとして準備できないからである。ここで、3.3 節で述べる混合音からのスペクトル包絡推定手法を使用する。そして、推定されたスペクトル包絡を用いて、変換元スペクトルを最もよく表現する SIR の値

変換先歌手の歌声テンプレート

を、W-PST 法により計算する（図 3 (d)）。

合成されたテンプレート

次に、その変換元歌手をモデル化した歌声スペクトルテンプレートを、変換先の歌手をモデル化した歌声スペクトルテンプレートに置き換える（図 3 (e)）。ここで変換先の歌手の歌
声スペクトルテンプレートは、単独歌唱の変換先音響信号から文献1) で提案された手法に(g) スペクトルを操作より推定する（図 3 (f)）。これにより、スペクトルの各周波数ビン（離散的に計算された周返還後のスペクトルテンプレート

波数成分）ごとに、変換元の歌手から変換先の歌手へ声質変換することで、どの程度パワーを変化させる必要があるかが計算できる。ここで、ノイズスペクトルテンプレートとしては共通のものを使用しているため、歌声が優勢でない周波数帯域はパワーを変化させる必要が返還後のスペクトル(h) 逆連続ウェーブ
レット変換
図3

無くなり、結果として伴奏音の音質は保存される。そして、変換元のスペクトルの各周波数
出力音響信号ビンのパワーを実際に変化させることで、変換先の歌声の声に変換された変換後のスペクト

声質変換処理の概要．

ルを得ることができる（図 3 (g)）。最後に，変換後のスペクトルに対して逆連続ウェーブ
レット変換（ICWT）をかけることで音響信号を再合成する（図 3 (h)）．このとき，位相

歌声とそれ以外の音に分けて考え，それぞれが別々の確率的スペクトルテンプレートからは元の観測スペクトルのものをそのまま利用する．

独立に生成され，足しあわされることで観測スペクトルが生成されたと考える．それらの 2

2.2 定 式 化

つのスペクトルテンプレートの加算の際に重みパラメータを導入し，重み付きで加算する

前述の手法の具体的な定式化は下記のようになる．

ことで，様々な SIR のスペクトルを表現できる．さらに，歌声スペクトルテンプレートは，

2.2.1 ウェーブレット変換によるスペクトルの計算

歌声包絡テンプレート（図 2 (a)）と駆動音源関数（図 2 (b)）の積によって生成されると

まず，入力音響信号に対して CWT をかけることで，スペクトログラムを計算する．本仮定する．この仮定はソースフィルタモデルを近似的に表現したものである．

研究では，マザーウェーブレットとしてガボールウェーブレットを用いる．ガボールウェー

図 3 に W-PST 法を用いた声質変換の概要を示す．まず入力として変換元音響信号と変ブレットによる CWT は，入力音響信号を x(t) とすると，下記のように定義される．

3c 2010 Information Processing Society of JapanVol.2010-MUS-86 No.7
2010/7/28

情報処理学会研究報告
IPSJ SIG Technical Report ∞1

W (b, a) = 
Ψ(t) = √

|a|

1
2πσ 2


x(t)Ψ

−∞



exp −

2

t
2σ 2





t−b
dt
a

(1)

exp (iω0 t)(2)

ここで，N (μ, σ 2 ) は，平均 μ，分散 σ 2 の正規分布を表す．さらに，調波構造を持つ歌声を表現する確率変数 Yv,f は，次式のように，スペクトル包絡の確率モデルと調波構造を表現するスペクトルの加算で表現できると仮定する（図 2）．3 節で述べたように，これはソースフィルタモデルの近似的表現である．

ただし，Ψ(·) は，Ψ(·) の共役複素数を表す現在の実装では、σ 2 を 3.75[ms] に、ω0 を 1 に
Yv,f = Yv,f
+ H(f ; f0 )設定している。b は時刻を表すパラメータで，W (b, a) は全ての b について（つまり，離散信(8)2∼ N (μv,f + H(f ; f0 ), σv,f)

Hz号の場合は全てのサンプルについて）計算される．a は周波数を表すパラメータで， 2πaw02exp(−(log f0 + log h − log f )2 /2θH)

H(f ; f0 ) = logに対応する．また，式 (1) の積分を無限の範囲で計算することは不可能なので，ウェーブ

(9)
(10)

h

レット Ψ(t) の区間を適当に切り詰めて計算する．本稿の実験では，−3σ < t < 3σ の範囲で計算をした．また，実装の際には，式 (1) は a の値を固定すると畳み込み演算であること2
ここで，Yv,f∼ N (μv,f , σv,f
) は歌声のスペクトル包絡を表現する確率変数であり，歌声包を利用して，高速フーリエ変換 (FFT) を用いた畳み込みの高速計算法を利用した．

絡テンプレートと呼ぶ．また，H(f ; f0 ) は F0 の値が f0 の声帯振動のスペクトルを表現し，次節以降で述べるテンプレートの推定処理では，計算時間を削減するために，10ms 間隔

駆動音源関数と呼ぶ（図 2 (b)）．なお，駆動音源関数 H(f ; f0 ) は確率変数ではないことにの離散的な b について（以降，フレームと呼ぶ）計算する．以降の処理は，それぞれの離散

注意が必要である．ただし，F0 と歌声包絡テンプレートとノイズスペクトルテンプレート

的な b の値について独立に行われるため，b の表記は省略し，対数パワースペクトル y(f ) を2
2，μn,f ，σn,fの推定方法は次節で述べるため、本節では既知のものパラメータ μv,f ，σv,fy(f ) = log (|W (b, a)|)と表記する．ただし，a と f には

f = logという関係がある．

2πa
w0

(3)

2
のとする。また、現在の実装では、θH
は、15 cent に設定している。

(4)ように表される．

以上をまとめると，歌声と伴奏音が混ざったスペクトルを表現する確率変数 Yf は下記の
Yf = log(exp(Yv,f+ H(f ; f0 ) + gv ) + exp(Yn,f + gn ))(11)2.2.2 確率的スペクトルテンプレート

確率変数 Yf はパラメータ (gv , gn ) に依存する．以降の説明では，便宜的に確率変数 Yf が歌声を含む混合音の対数パワースペクトル y(f ) は，ある確率変数（の集合）Yf から生従う確率密度関数を pf (y; gv , gn ) と記す．

成されると仮定する．この確率変数 Yf を確率的スペクトルテンプレートと呼ぶ．次に，Yf

2.2.3 スペクトルテンプレートの加算の近似は次式により 2 つの異なるスペクトルテンプレート Yv,f と Yn,f に分割できると仮定する．

式 (11) で表される確率的スペクトルテンプレート Yf の確率密度関数は，解析的に計算す(5)

ることは困難であるので，正規分布を用いて近似計算する．関数 l(x1 , x2 ) = log(exp(x1 ) +ただし，Yv,f は歌声のスペクトルを表し，歌声スペクトルテンプレートと呼ばれ，Yn,f はYf = log(exp(Yv,f + gv ) + exp(Yn,f + gn ))exp(x2 )) の (x1 , x2 ) = (μv,f + H(f ; f0 ) + gv , μn,f + gn ) における 1 次のテーラー展開は

歌声以外の音（伴奏音）のスペクトルを表し，ノイズスペクトルテンプレートと呼ばれる．

l(x1 , x2 ) ≈gv と gn はそれぞれのテンプレートの重みであり，それらを変化させることで歌声とその他の音の SIR を変化させることができる．なお，式 (5) においては，パワースペクトルの加

exp (μv,f + H(f ; f0 ) + gv )x1
exp (μv,f + H(f ; f0 ) + gv ) + exp (μn,f + gn )+法性を仮定している．

Yv,f と Yn,f が，次式のように，（対数周波数軸上で）正規分布に従うと仮定する．

exp (μn,f + gn )x2 + Cexp (μv,f + H(f ; f0 ) + gv ) + exp (μn,f + gn )(12)のように計算される．ただし，C は x1 と x2 とは独立な定数である．ここで，パラメータ2
Yv,f ∼ N (μv,f , σv,f)

(6)

gv ，gn が固定された場合，式 (12) が x1 と x2 の重み付き加算であることに注意すると，

2
Yn,f ∼ N (μn,f , σn,f
)

(7)


Yf = l(Yv,f+ H(f ; f0 ) + gv , Yn,f + gn ) が従う確率密度関数 pf (y; gv , gn ) は，4

c 2010 Information Processing Society of JapanVol.2010-MUS-86 No.7
2010/7/28

情報処理学会研究報告
IPSJ SIG Technical Reportpf (y; gv , gn ) ≈ N (y; μf (θv , θn ), σf2 (θv , θn ))μf (gv , gn ) = log(exp (μv,f + H(f ; f0 ) + gv ) + exp (μn,f + gn ))22+ (exp (μn,f + gn ))2 σn,f(exp (μv,f + H(f ; f0 ) + gv ))2 σv,fσf2 (gv , gn ) =(exp (μv,f + H(f ; f0 ) + gv ) + exp (μn,f + gn ))2(13)

る．時刻 b におけるフィルターを ζ(b, a) と書くと，スペクトルの位相は元のものを使うの(14)で，変換後のウェーブレットスペクトログラム Ŵ (b, a) はW (b, a)Ŵ (b, a) =(|W (b, a)| + ζ(b, a))(19)|W (b, a)|で表される．ただし，a と f には，式 (4) のような関係がある．ウェーブレットスペクトロ(15)

グラムを時間信号 x̂(t) に変換する ICWT は，次式で定義される．

のように表現される．ただし，N (y; μ, σ 2 ) は，平均 μ，分散 σ 2 の正規分布の確率密度関数x̂(f ) =

を表す．

1
CΨ

 ∞ ∞−∞1t−bŴ (b, a)  Ψ
a
|a|
−∞



1
dadb
a2

(20)

2.2.4 準ニュートン法によるパラメータ最適化

ただし，CΨ は定数であるが，全ての時刻で同じ値をとるため厳密に計算する必要はない．

SIR を表すパラメータ (gv , gn ) の最適化には，BFGS（Broyden-Fletcher-Goldfarb-

ICWT の計算においても CWT と同様で，FFT による畳み込み演算の高速計算法を用いてShanno）公式に基づく準ニュートン法を使用する．準ニュートン法は山登り法の一種で実装した．

あり，反復的にパラメータを更新する．本モデルにおいて，最小化すべき目的関数 Q(gv , gn )3. 歌声包絡テンプレートの推定は，log N (y(f ); uf (gv , gn ), σf2 (gv , gn ))Q(gv , gn ) = −(16)

前章では，変換前の歌声の歌声包絡テンプレートとノイズスペクトルテンプレート，およfび変換先の歌声の歌声包絡テンプレートが与えられているという条件で，歌声変換手法につで表される．ただし，y(f ) は観測スペクトルである．

いて議論した．本節では，それらのテンプレートの具体的な構成方法と，テンプレートを入

2.2.5 ウェーブレット変換に基づく声質の変換

力音響信号から推定する手法を述べる．

以上により，y(f ) を最もよく表現する重み gv と gn の値が推定でき，その時の合成後の3.1 混合回帰モデルによるテンプレートの表現

スペクトルテンプレートの確率密度関数 pf (y; gv , gn ) が計算できる．次に，変換元の歌声包スペクトルテンプレートを表現するモデルとして，文献 1) と同様に，各回帰要素として

22を変換先の歌声包絡テンプレート μ̂v,f と σ̂v,fと絡テンプレートのパラメータ μv,f と σv,f

線形回帰を使用した混合回帰モデル14) を導入する．前章で述べたように，本手法において置き換えて，変換先のスペクトルテンプレート pˆf (y; gv , gn ) を計算する．スペクトル y(f )ˆ ) へ，下記の式により変換する．
を新しいスペクトル y(fŷ(f ) = y(f ) + ζ(f )(17)

ζ(f ) = Ey [p̂f (y; gv , gn )] − Ey [pf (y; gv , gn )](18)はスペクトルテンプレートはある周波数 f における対数パワーの分布が正規分布で表現されるモデルを用いて定義される必要があるが，このモデルはその用件を満たしている．混合回帰モデルは任意の非線形回帰を複数の線形回帰によって近似するモデルで，スペクトル包絡の形状について仮定を置かず，学習データのみに基づいてスペクトル包絡を推定する．混

ただし，E[·] は期待値を表す．ζ(f ) はフィルターの役割を果たす関数で，元のスペクトルを合数 M の混合回帰モデルは，m(1, · · · , M ) を各線形回帰モデルのインデックスとすると、変換先の歌手の歌声に変換するために操作が必要な周波数帯域とその操作量を表している．

パラメータとして、それぞれの線形回帰モデルの傾き av,m と切片 bv,m 、各線形回帰モデル．以
また，歌声の音量を調整したい場合は，パラメータ ĝv を増減させることで実現できる．

2をとる。ゲート関数としの守備範囲を決めるゲート関数のパラメータ ψv,m 、μv,m 、σv,m上により変換後のスペクトル ŷ(f ) を得ることができる．

ては、次式で定義される正規化ガウス関数15) を用いた．
2
ψv,m N (f ; μv,m , σv,m)2
Gm (f ; ψv,m , μv,m , σv,m)=M2ψ  N (f ; μv,m , σv,m)m =1 m最後に，得られたスペクトルを逆ウェーブレット変換して，変換後の音響信号を得る．前
述のように，計算時間の削減のため，式 (17) の計算は 10ms のフレーム毎に行われるので，上述の ŷ(f ) はその他の b の値では計算されない．そこで，隣り合うフレーム間の ζ(f ) のここで、ψv,m は各ガウス関数の重みを決めるパラメータで、ψv,m ≥ 0 かつ

(21)
M
ψ
=1
m=1 v,m2
である．また、μv,m と σv,mは、ガウス関数の平均と分散である。このモデルでは、歌声値を線形補間することで，全ての b について ζ(f ) を計算し，式 (17) により ŷ(f ) を計算す

5c 2010 Information Processing Society of JapanVol.2010-MUS-86 No.7
2010/7/28

情報処理学会研究報告
IPSJ SIG Technical Report2
包絡テンプレートのパラメータである平均 μv,f と分散 σv,fはStep 1 混合回帰モデルのパラメータを EM 法により推定する．
Step 2 ki を次式により更新する．
Hi μv,fi,h − yi,h
2
h=1
σv,f
i,h
ki =
1
Hi
h=1 σ 2v,fi,hM

μv,f =

2
Gm (f ; ψv,m , μv,m , σv,m)(av,m f + bv,m )

(22)

2
2
Gm (f ; ψv,m , μv,m , σv,m)2 βv,m

(23)

m=1
M
2
σv,f
=
m=1

(26)

として表現する．ただし，M は混合数を表す．現在の実装では、M を 10 に設定している。

Step 3 1 に戻る．

このモデルの未知パラメータは、EM（Expectation and Maximization）法により推定すki 以外のパラメータの初期値として，周波数軸の定義域 (現在の実装では 60Hz〜7500Hz)

ることが可能である．ノイズスペクトルテンプレートについても同様で，未知パラメータをを M 等分し，m 番目の分割について，(fi,h , yi,h ) の回帰係数を計算したものを am と bm

2
2
{ψn,m , μn,m , σn,m, an,m , bn,m , βn,m
} と置き，同様の形式で表現する．

2
1の初期値に，fi,h の平均と分散を μm と σmの初期値に設定し，ψm の初期値は Mとした．

3.2 単独歌唱からのテンプレート推定

ノイズスペクトルテンプレートについては，si (i = 1, . . . , I) を調波構造でなくスペクトル
そのものと考えることで，同様に推定できる．現在の実装では、Step 1 の EM 法の反復回

単独歌唱の音響信号が与えられている場合は，歌声包絡テンプレートとノイズスペクトル数は 1 にし、Step 0〜3 全体の反復回数を 30 回に設定している。

テンプレートは，個別に推定する．歌声包絡テンプレートは、各母音毎に独立に推定され、例えば母音/a/のテンプレートを推定する際は、学習データ中の/a/のラベルが付与されて3.3 混合音からのテンプレート推定

いるフレームのみを用いて推定する。ノイズスペクトルテンプレートは全体で 1 つが推定

混 合 音 か ら テ ン プ レ ー ト を 推 定 す る 場 合 は ，歌 声 包 絡 テ ン プ レ ー ト と ノ イ ズ

される。現在の実装では、ノイズスペクトルテンプレートの推定には、歌声を含まない伴奏

ス ペ ク ト ル テ ン プ レ ー ト を 同 時 に 推 定 す る 必 要 が あ る ．I 個 の 観 測 ス ペ ク ト ル

のみの音響信号（カラオケトラック）を使用している。

y1 (f ), · · · , yi (f ), · · · , yI (f ) を観測したとする．推定すべき歌声テンプレートのパラメー

複数の調波構造からその元となるスペクトル包絡を推定する場合，フレームごとの音量の

2
2, av,m , bv,m , βv,m
} とし，ノイズテンプレートのパラメータはタは θv = {ψv,m , μv,m , σv,m違いを考慮に入れる必要がある．そのため，本研究では各フレームの音量を正規化するため

22
θn = {ψn,m , μn,m , σn,m, an,m , bn,m , βn,m} とする．i 番目のスペクトルにおける駆動音源のパラメータを導入し，それも未知パラメータとして推定することでこの問題を解決する．

関数を加えた後の歌声スペクトルテンプレートは，

μv,f,i = μv,f + H(f ; f0 (i))

学習データとして与えられた I フレーム分の調波構造 si (i = 1, . . . , I) の h 次倍音の周波数 fi,h とその対数パワー yi,h が，

si = {(fi,1 , yi,1 ), . . . (fi,h , yi,h ), . . . (fi,Hi , yi,Hi )}(24)

いて既知であるとする．
前章では，対数正規分布の加算を 1 次のテイラー展開を用いて近似計算した．しかし，として表されるとする．この時，最大化したい尤度関数は，次式で表される．
I

i

Hi

2
log N (yi,h + ki ; μv,fi,h , σv,f
)
i,h

(27)

と表すことができる．ただし，i 番目の観測スペクトルの F0 である f0 (i) は全ての i につ得られた式 (13)〜(15) は複雑な形状となり，θv ，θn を最適化するのは困難である．そこ(25)

で本節では，対数正規分布の加算を定義に従って厳密に計算した後，パラメータを近似的

h

ここで，ki は各調波構造の音量をフレーム間で正規化するオフセットパラメータである．混に推定するというアプローチをとる．合成後のスペクトルテンプレートの確率密度関数を

合回帰モデルのパラメータと ki を同時に最適化することは困難なので，それらを反復的に

pi,f (y; θv , θn , gi,v , gi,n ) 1 と書くと，目的関数 L は，

更新していく．
パラメータは下記の手続きで推定される．

1 2.2.3 節と異なり，観測するスペクトルの番号 i ごとに確率密度関数の形状が異なるので，添字 i を追加している．

Step 0 ki = 0 とし，その他のパラメータに対して後述のように初期値を与える．

6c 2010 Information Processing Society of JapanVol.2010-MUS-86 No.7
2010/7/28

情報処理学会研究報告
IPSJ SIG Technical Report

すると，πi,r,f と log(exp(yi (f )) − exp(Ui,r,f )) は定数となるため，式 (31) を用いて，gi,vIL=log pi,f (y; θv , θn , gi,v , gi,n )df

i=1

 yi (f )

I

=

log(
−∞

i=1

 yi (f )

I

log(

=

−∞

i=1

と θv を最適化できる．また，gi,v と θv を固定した場合も同様で，式 (29) からサンプリン
グにより式 (31) と同様の式を導出し，gi,n と θn を最適化する．
しかし，式 (31) は和の対数の形をしているため，未だ直接の最適化が困難である．そこ2
N (log(exp(yi (f )) − exp(U )); μv,f,i + gi,v , σv,f)で，EM アルゴリズムに似た反復法によって，式 (31) を反復的に最適化する．便宜的に，推2)
N (U ; μn,f + gi,n , σn,f



(28)

exp(yi (f ))dU )dfexp(yi (f )) − exp(U )定したいパラメータを λ = {gi,v , θv } と書く．また，一回前の反復におけるパラメータの推(29)

定値を λ と置く．まず，変数 zi,r,f を導入し，

2
N (log(exp(yi (f )) − exp(U )); μn,f + gi,n , σn,f)2)
N (U ; μv,f,i + gi,v , σv,fexp(yi (f ))dU )dfexp(yi (f )) − exp(U )2)πi,r,f N (log(exp(yi (f )) − exp(Ui,r,f )); μv,f,i + gi,v , σv,f
R
2


π
N
(log(exp(y
(f
))
−
exp(U
));
μ
+
g
,
σv,f
)
i
i,v
v,f,i
i,r ,fr  =1 i,r ,fzi,r,f =

(33)

(30)


と置き，λ を用いた計算した zi,r,f を zi,r,fと書く．このとき，zi,r,f を固定し新たな目的

関数 Q1 (λ|λ )と表現される．ここで，gi,v と gi,n は，3.2 節の ki と同様で，音量をフレーム間で正規化

するオフセットパラメータである．また，本節では，歌声包絡テンプレートとノイズスペ


I

R

zi,r,flog πi,r,fQ1 (λ|λ ) =クトルテンプレートの SIR を調整する役割も持っている．実際の実装では，連続ウェーブ

i=1 r=1

レット変換は周波数軸に対して離散的に計算しているため，f に関する積分は和の演算で置

2
N (log(exp(yi (f )) − exp(Ui,r,f )); μv,f,i + gi,v , σv,f)df

き換えられる．
ここで推定すべきパラメータは {gi,v , gi,n , θv , θn } である．これらのパラメータを全て同(34)

を λ に関して最適化する操作と，最適化された λ を用いて zi,r,f を再計算する操作を反復時に最適化するのは困難であるので，逐次的に最適化する．まず，gi,n と θn を固定して，すると真の目的関数 L が最大化できる．証明は付録を参照されたい．

式 (29) による gi,v と θv の最適化と，gi,v と θv を固定して，式 (30) による gi,n と θn の最式 (34) をよく見ると，πi,r,f は最適化に無関係であることがわかり，関数 Q2 (λ|λ )適化を交互に繰り返すことを考える．まず，gi,n と θn を固定して考えると，式 (29) の和のQ2 (λ|λ ) =内部は期待値の計算と考えることができる．そこで，U を期待値の計算をサンプリングに



I

R


2
zi,r,flog N (log(exp(yi (f )) − exp(Ui,r,f )); μv,f,i + gi,v , σv,f)df
i=1 r=1

より和の計算で近似することにより，gi,v と θv の近似的な最適化を可能にする．具体的に

(35)



2
) を U = yi (f ) で切断した，上に有界な単一切断正規は，正規分布 N (U ; μn,f + gi,n , σn,fの最適化は，Q1 (λ|λ ) の最適化と等価であることがわかる．さらに，Q2 は定数項 z の存在

分布からそれぞれの i，f について R 個ずつのサンプル (Ui,1,f , · · · , Ui,r,f , · · · , Ui,R,f ) をを除くと，式 (25) と同様の形式をしていることがわかる．すなわち，Q2 は 3.2 節で述べた

サンプリングしたとき，目的関数 L は，単独歌唱からのテンプレート推定の場合と同様に最適化できることがわかる．



I

L≈
πi,r,f =Step 0 gi,v = 0，gi,n = 0 とし，その他のパラメータに対して後述のように初期値を与2πi,r,f N (log(exp(yi (f )) − exp(Ui,r,f )); μv,f,i + gi,v , σv,f
) (31)

log
i=1以上をまとめるとパラメータは下記の手続きで推定される．

R

える．

r=1

exp(yi (f ))(exp(yi (f )) − exp(Ui,r,f )) R

yi (f )−∞

2
N (U ; μn,f + gi,n , σn,f)dU

Step 1 gi,n と θn を固定して，式 (29) の U をサンプリングする．

(32)

Step 2 サンプリングした U と現在のパラメータ gi,v ，θv を用いて，式 (33) の zi,r,f を計算する．

と近似できる．現在の実装では R の値を 300 に設定している。ここで，gi,n と θn を固定

7c 2010 Information Processing Society of JapanVol.2010-MUS-86 No.7
2010/7/28

情報処理学会研究報告
IPSJ SIG Technical Report

Step 4 Step 2〜3 の反復が規定回数を超えた場合は Step 5 へ，そうでない場合は Step 2

Log Power

Log Power

節の反復的な最適化法を利用する．

Log Power

Step 3 計算された zi,r,f を用いて，式 (35) の Q2 関数を最適化する．この最適化には 3.2に戻る．

Step 5 gi,v と θv を固定して，式 (30) の U をサンプリングする．

500

Step 6 サンプリングした U と現在のパラメータ gi,n ，θn を用いて，式 (33) の zi,r,f を

2000
1000
Frequency (Hz)4000 7000(a) 3.2節の手法による単独歌唱からの
推定結果

計算する．

Step 7 計算された zi,r,f を用いて，式 (35) の Q2 関数を最適化する．この最適化には 3.2

図4

節の反復的な最適化法を利用する．

500

2000
1000
Frequency (Hz)4000 7000(b) 提案法による推定結果

500

2000
1000
Frequency (Hz)4000 7000(c) 混合音から抽出した調波構造からの3.2節の手法による推定結果

混合音からのスペクトル包絡推定の例．
「RWC 研究用音楽データベース：ポピュラー音楽（RWC-MDB-P2001）」16) の No.7 の楽曲から推定した/i/の音素のスペクトル包絡である．

Step 8 Step 2〜3 の反復が規定回数を超えた場合は Step 9 へ，そうでない場合は Step 6に戻る．

でいることが見てとれる。これにより、提案法が伴奏音の影響を低減できることがわかる。

Step 9 Step 1〜8 の反復が規定回数を超えた場合は終了する．そうでない場合は Step 1

次に，声質変換によるスペクトル変化の実例を示す．図 5 は，No.7 の楽曲（図 5(a)）の声質を，ボーカルをキャンセルした場合（図 5(b)），No.13 の歌手の声に変換した場合（図に戻る．
歌声包絡テンプレートの初期値は、今回の推定対象とは異なる歌手の単独歌唱の音響信号か

5(c)），No.20 の歌手の声に変換した場合（図 5(d)）のスペクトルの変化である．また，図ら、ノイズスペクトルテンプレートの初期値は，歌声の入っていない音楽音響信号（カラオ

6 は，図 5 において元の楽曲を No.20 の声質に変換した場合の，それぞれに対応するスペケトラック）から、それぞれ 3.2 節の手法により推定したパラメータの値を使用する．

クトルテンプレートの例である．なお、これらの使用した楽曲は全て女性で、異なる歌手

4. 実のものである。また、1 楽曲あたり各母音が 2000〜5000 フレーム程度含まれている。図中

装のスペクトルには音素/i/の音が含まれている．なお，ボーカルキャンセルとは，2.2.5 節の

上記の技術を用い，混合音中の歌声の声質変換を実装した．声質の変換は正解が存在しな

ĝv を −∞ に設定して声質を変換した場合であり，声質を変換するのではなく歌声の音量をい操作であり，定量的な評価が困難であるので，ここでは声質を変換した場合の実験結果例

下げる変換に相当する．図より，伴奏音に起因する周波数成分は変化していないが，400Hz

をいくつか紹介する．被験者実験等による評価実験を行うことは今後の課題となる．

付近のピークや，2500〜4500Hz 付近のピークなど，歌声の周波数成分の形状が変化してい「RWC 研究用音楽
まず，混合音からのスペクトル包絡の推定の実行例を示す．図 4 は，
16)

データベース：ポピュラー音楽（RWC-MDB-P-2001）」

ることがわかる．特にボーカルをキャンセルした場合は，2500〜4500Hz 付近のピークが顕の No.7 の楽曲について，歌声

著に無くなっている．

のスペクトル包絡を，単独歌唱から 3.2 節の手法を用いた推定したもの（図 4(a)），混合音

ここで図示した以外にも，いくつかの歌手の組に対して変換を実行した．聴感上，ボーカ

から 3.3 節の手法を用いて推定したもの（図 4(b)），混合音から抽出した調波構造を用いてルキャンセルに関しては，わずかに歌声が残っているものの，伴奏音の音質には影響を与え

3.3 節の手法により推定したもの（図 4(c)）を図示している．(a) は単独歌唱から推定してずに，歌声の音量を低減できていた．声質変換に関しては，主観的な印象では，変換後もわいるので理想的な推定結果と考えることができ、提案法の推定結果 (b) がどれだけ (a) に近ずかに元の歌手の特徴が残りながらも，変換先の歌手の特徴を持った声に変換されている

いかが問題となる。(c) は、伴奏音の影響を考慮せず、伴奏音が重畳した状態から推定したように聞こえた．一方で，楽曲によっては，変換元の歌手の声と変換先の歌手の声が混ざっ場合である。図 4 から見てとれるように、(b) では全体に分散が大きくなる傾向や広域のパたような声になる場合もあった．また，異なる性別の歌手の声に変換する場合や，歌声の音

ワーの弱い部分で歪みが増える傾向はあるものの、(a) に近い推定結果が得られていること

量を大きく増加させた場合に，不自然な音声になることがあった．これは，元のスペクトルがわかる。一方、(c) では、伴奏音の影響により、(a) と比較してスペクトルが大きく歪んで伴奏音に埋もれてしまっている周波数帯域を無理に増大させたことにより，位相が不自然

8

c 2010 Information Processing Society of JapanVol.2010-MUS-86 No.7
2010/7/28

500

1000

2000

500

4000 7000

1000

2000

Frequency (Hz)Frequency (Hz)(a) 元のスペクトル（No.7）(b) ボーカルキャンセル

4000 7000

500

2000
1000
Frequency (Hz)Log Power

500

4000 7000

2000
1000
Frequency (Hz)4000 7000(b) No. 20の歌手の歌声包絡テンプレートを用いて操作されたスペクトルテンプレート(a) 元のスペクトル（No.7）を表現する合成されたスペクトルテンプレート
図6

Log Power

Log Power

Log Power

Log Power

Log Power

情報処理学会研究報告
IPSJ SIG Technical Report

声質変換の際のテンプレートの置換の例．
「RWC 研究用音楽データベース：ポピュラー音楽（RWC-MDB-P2001）」16) の No.7 の楽曲のスペクトルを表現する合成されたテンプレート (a) と，その歌声を No.20 のテンプレートに置換したテンプレート (b) を図示する．(a) は図 5 (a) に，(b) は図 5 (d) に対応する．

した状態で提供される一般の音楽音響信号に対して適用可能にした．本技術を実装し，実際にポピュラー音楽に対して実行することで，提案法により正しく声質が変換されることを確認した．今後の課題は，被験者を用いた評価実験を行い，提案法の性能を評価することであ

500

2000
1000
Frequency (Hz)(c) No. 13 の歌手に変換
図54000 7000

500

2000
1000
Frequency (Hz)4000 7000

る．また，本稿では歌声の音素と F0 のラベルが付与されていることを仮定し，母音に対してのみ処理をすることで、声質が変換できることを確認した。今後のさらなる性能向上のた(d) No. 20 の歌手に変換

めには、その仮定をなくし，全ての音素に対して処理をするために本手法を拡張していく予

声質変換によるスペクトル変化の例．
「RWC 研究用音楽データベース：ポピュラー音楽（RWC-MDB-P2001）」16) の No.7 の楽曲に対して，(b) ボーカルキャンセル，(c)No.13 の歌手の声に変換，(d)No.20 の歌手の声に変換の 3 種類の処理をした場合のスペクトルを図示する．図中の点線はスペクトルに含まれる基本周波数（約 490Hz）とその倍音周波数を表している。

定である。
謝辞 本研究の一部は CrestMuse プロジェクト (JST CREST) の支援を受けた．また、本稿での実験に「RWC 研究用音楽データベース：ポピュラー音楽（RWC-MDB-P-2001）」16)を使用した。

になったためだと考えられる．これに対しては，そのような周波数帯域では，歌声の周波数

参

成分を正弦波重畳モデル等で別に再合成して足しあわせるなどの処理が有効であると考え

考

文

献

1) 藤原弘将，後藤真孝，奥乃 博：多重奏中の歌声の基本周波数と音素を同時に推定可能な新たなフレームワーク，情報処理学会研究報告，Vol.2009-MUS-81 (2009).
2) Goto, M.: Active Music Listening Interfaces Based on Signal Processing, Proceedings of the 2007 IEEE International Conference on Acoustics, Speech, and Signal
Processing (ICASSP 2007), pp.IV–1441–1444 (2007).
3) Yoshii, K., Goto, M., Komatani, K., Ogata, T. and Okuno, H.G.: Drumix: An Audio Player with Real-time Drum-part Rearrangement Functions for Active Musicられる．

5. ま と め
本稿では，混合音中の歌声の声質変換を実現する手法について述べた．具体的には，W-

PST 法1) を応用して，歌声のみの周波数成分のみを操作することを可能にした．さらに，混合音の音響信号中の歌声のスペクトル包絡推定手法を開発することで，歌声と伴奏と混在

9c 2010 Information Processing Society of JapanVol.2010-MUS-86 No.7
2010/7/28

情報処理学会研究報告
IPSJ SIG Technical Report

Listening, IPSJ Journal, Vol.48, No.3, pp.1229–1239 (2007).
4) Itoyama, K., Goto, M., Komatani, K., Ogata, T. and Okuno, H.G.: Instrument
Equalizer for Query-by-Example Retrieval: Improving Sound Source Separation
based on Integrated Harmonic and Inharmonic Models, Proceedings of the 9th International Conference on Music Information Retrieval (ISMIR 2008), pp.133–138(2008).
5) Bonada, J., Celma, O., Loscos, A., Ortola, J., Serra, X., Yoshioka, Y., Kayama, H.,
Hisaminato, Y. and Kenmochi, H.: Singing Voice Synthesis Combining Excitation
plus Resonance and Sinusoidal plus Residual Models, Proceedings of International
Computer Music Conference (2001).
6) 剣持秀紀，大下隼人：歌声合成システム VOCALOID，情報処理学会研究報告，
Vol.2007-MUS-72, pp.25–28 (2007).
7) 酒向慎司，宮島千代美，徳田恵一，北村 正：隠れマルコフモデルに基づいた歌声合成システム，情報処理学会論文誌，Vol.45, No.3, pp.719–727 (2004).
8) Sinsy - HMM-based Singing Voice Synthesis System: http://www.sinsy.jp/.
9) Stylianou, Y., Cappé, O. and Moulines, E.: Continuous probabilistic transform forvoice conversion, IEEE Transactions on Speech and Audio Processing, No.2, pp.
131–142 (1998).
10) Mouchtaris, A., der Spiegel, J.V. and Mueller, P.: Nonparallel training for voice
conversion based on a parameter adaptation approach, IEEE Transactions on Audio, Speech and Language Processing, Vol.14, No.3, pp.952–963 (2006).
11) Toda, T., Black, A.W. and Tokuda, K.: Voice conversion based on maximum likelihood estimation of spectral parameter trajectory, IEEE Transactions on Audio,
Speech and Language Processing, Vol.15, No.8, pp.2222–2235 (2007).
12) 河原英紀，生駒太一，森勢将雅，高橋 徹，豊田健一，片寄晴弘：モーフィングに基づ
く歌唱デザインインタフェースの提案と初期検討，情報処理学会論文誌，Vol.48, No.12,
pp.3637–3648 (2007).
13) 大西壮登，高橋 徹，入野俊夫，河原英紀：一般逆行列を用いた母音情報に基づく声質変換法について，電子情報通信学会技術報告，No.282, pp.75–80 (2007).
14) Jacobs, R.J., Jordan, M., Nowlan, S.J. and Hinton, G.E.: Adaptive mixtures of
local experts, Neural Computation, Vol.3, pp.79–87 (1991).
15) Xu, L., Jordan, M. I. and Hinton, G. E.: An alternative model for mixtures ofexperts, Advances in Neural Information Processing Systems 7, pp.633–640 (1994).
16) 後藤真孝，橋口博樹，西村拓一，岡 隆一：RWC 研究用音楽データベース:研究目的で利用可能な著作権処理済み楽曲・楽器音データベース，情報処理学会論文誌，Vol.45,
No.3, pp.728–738 (2004).

付録 3.3 節の反復アルゴリズムの妥当性の証明
証明. 3.3 節において，式 (31) の L は，Jensen の不等式より



 R

I

L(λ) =


zi,r,f

log



i=1
I

r=1
R

zi,r,flog≥
i=1 r=12πi,r,f N (xi,r,f ; μv,f,i + gi,v , σv,f)
zi,r,f

df

(36)

2
)
πi,r,f N (xi,r,f ; μv,f,i + gi,v , σv,fdf = F (λ|λ )zi,r,f(37)と変形できる．ただし，xi,r,f は，

xi,r,f = log(exp(yi (f )) − exp(Ui,r,f ))である．このとき，L(λ) − L(λ ) = F (λ|λ ) − F (λ |λ ) +



(38)



I

zi,r,flog
i=1zi,r,f
zi,r,f
df(39)

が成立する．右辺第三項は非負なので F (λ|λ(−1) ) の λ に関する最大化は目的関数 L(λ) を増加させることがわかる．さらに，






I





zi,r,flog zi,r,fdf

Q1 (λ|λ ) = F (λ|λ ) +

(40)

i=1

と変形でき，右辺第二項は λ に無関係な項であるので，F (λ|λ ) の λ に関する最大化は，

Q1 (λ|λ ) の λ に関する最大化と等価である．以上より，Q(λ|λ ) を最大化させることで，目的関数 L(λ) が増加することが示される．

10c 2010 Information Processing Society of Japan