INTERSPEECH 2013

Probabilistic speech F0 contour model incorporating statistical vocabularymodel of phrase-accent command sequence
Tatsuma Ishihara† ,Hirokazu Kameoka†‡ , Kota Yoshizato† , Daisuke Saito† , Shigeki Sagayama† ,†Graduate School of Information Science and Technology, The University of Tokyo, Japan‡NTT Communication Science Laboratories, NTT Corporation, Japan{ishihara,yoshizato,kameoka,dsaito,sagayama}@hil.t.u-tokyo.ac.jp

Abstract

ponents, one important challenge is to solve an inverse probelmof estimating the Fujisaki-model parameters automatically froma raw F0 contour.
However, this problem has been a difﬁcult task. Several
techniques have already been developed [2, 3, 4], but so far
with limited success due to the ill-posedness of the inverse
problem and the analytical complexity of the Fujisaki model.
We have previously derived a stochastic model of speech F0
contours by translating the Fujisaki model into a probabilistic
generative model [5, 6, 7]. This reformulation has successfully allowed us to derive an efﬁcient algorithm based on the
expectation-maximization (EM) algorithm for estimating the
Fujisaki-model parameters from a raw F0 contour.
A key idea in our previous model is that the sequence ofthe phrase and accent command pair (i.e., the parameters that
we want to estimate) is modeled as a path-restricted hidden
Markov model (HMM) so that estimating the state transition of
the HMM directly amounts to estimating the Fujisaki-model parameters. Generally speaking, the better top-down knowledge
is incorporated, the better the solution to an inverse problem isobtained. In speech recognition, for example, designing an appropriate state transition network of an HMM allows us to effectively search for a linguistically likely phoneme sequence thatbest explains an observed sequence. With the same strategy, we
would want to design an appropriate state transition network of
our HMM to effectively reduce the solution space of the phrase
and accent command sequences. Here, a question is how to design a state transition network for our model. In normal speech,
many phrases or sentences share the same intonation pattern.
This is because an intonation pattern is usually determined by
the grammatical structure of an uttered sentence or the accent
type each phrase is associated with. Thus, it may be natural
to hypothesize that phrase and accent command sequences aregoverned by a vocabulary model. More speciﬁcally, we assume
that we have a dictionary consisting of a ﬁnite number of leftto-right HMM templates and a sequence of the phrase and accent command pair is generated according to a concatenationof those templates. In this paper, we formulate a probabilistic
F0 contour model by designing a state transition network of our
HMM based on the above assumption and derive a parameter
optimization algorithm.

We have previously proposed a generative model of speech F0
contours, based on the discrete-time version of the Fujisaki
model (a model of the mechanisim for controlling F0 s throughlaryngeal muscles). One advantage of this model is that it allows us to apply statistical methods to estimate the Fujisakimodel parameters from speech F0 contours. This paper proposes a new generative model of speech F0 contours incorporating a vocabulary model of intonation patterns. A parameter
inference algorithm for the present model is derived. We quantitatively evaluated the performance of our parameter inference
algorithm.

1. Introduction
All human speech has expression. We recognize it as part of
the humanness of speech, and it is a quality listeners expect toﬁnd in daily human communication. Without expression speechsounds lifeless and artiﬁcial. Our future goal is to quantify theexpression of speech in such a way that we can analyze andsynthesize it based on an engineering approach.
Fundamental frequency (F0 ) contours in normal speech
convey various types of non-linguistic information such as
speaker’s identity, emotion, attitude and intention. Modelingthe F0 contours of speech utterances can thus be potentially
useful for many speech applications, including emotion recognition, speaker recognition, speech synthesis, and dialogue systems, where prosodic features play an important role. F0 contours consist of two major components: major-scale pitch variations over the duration of the prosodic units (called the phrase
component) and smaller-scale pitch variations in accented syllables (called the accent component). This is justiﬁed by thefact that the thyroid cartilage involves two mutually independent types of movement with different muscular reaction times.
Speciﬁcally, the phrase and accent components respectively
correspond to contributions associated with the translation androtation movements of the thyroid cartilage. The Fujisaki model[1] is a well-founded mathematical model that describes an F0
contour as the sum of these two contributions. This model isknown to approximate actual F0 contours of speech surprisingly well when the model parameters are chosen appropriately,
and its validity has been shown for many, typologically diverse
languages. For this reason, and thanks to the intuitive association of the model parameters with the mechanical factors in thecontrol mechanism of phonation, the Fujisaki model has been
widely used with notable success to design F0 contours for synthesizing natural speech.
Since a prosodic feature in speech is predominantly characterized by the levels and timings of the phrase and accent com-

Copyright © 2013 ISCA

2. Probabilistic Pitch Contour Model
2.1. Original Fujisaki Model
The Fujisaki model [1] assumes that an F0 contour on a logarithmic scale, x(t), where t is time, is the superposition of threecomponents: a phrase component xp (t), an accent component

1017

25- 29 August 2013, Lyon, France

Figure 1: Original Fujisaki model [1].

xa (t), and a base component xb :x(t) =xp (t) + xa (t) + xb .

(1)
Figure 2: Previous HMM topology for phrase-accent command se-

The phrase component xp (t) consists of the major-scale pitch
variations over the duration of the prosodic units, and the accent
component xa (t) consists of the smaller-scale pitch variationsin accented syllables. These two components are modeled as
the outputs of second-order critically damped ﬁlters, one being excited with a command function up (t) consisting of Diracdeltas (phrase commands), and the other with ua (t) consistingof rectangular pulses (accent commands):
xp (t) =Gp (t) ∗ up (t),(α2 te−αt
Gp (t) =0xa (t) =Ga (t) ∗ ua (t),(β 2 te−βt
Ga (t) =0

quence modeling [5, 6, 7]. In state r0 , µp [k] and µa [k] are both constrained to be zero. In state p1 , µp [k] can take a non-zero value, Ap [k],whereas µa [k] is still restricted to zero. In state p1 , no self-transitions
are allowed. In state r1 , µp [k] and µa [k] again become zero. This
path constraint restricts µp [k] to consisting of isolated deltas. State a0
leads to states a1 , . . . , aN , in each of which µa [k] can take a different(n)non-zero value Aa , whereas µp [k] is forced to be zero. Direct statetransitions from state an to state an′ without passing through state r1
are not allowed. This constraint restricts µa [k] to consisting of rectangular pulses.

(2)(t ≥ 0),(t < 0)(3)
Output sequence: {o[k]}K
k=1
State sequence: {sk }K
k=1
Output distribution: P (o[k]|sk ) = N (o[k]; csk , Υsk )Mean sequence: —[k] = (µp [k], µa [k])T = csk
Transition probability: ϕi′ ,i = log P (sk = i|sk−1 = i′ )(4)(t ≥ 0),(t < 0)(5)where ∗ denotes convolution over time. The baseline component xb is a constant value related to the lower bound of the
speaker’s F0 , below which no regular vocal fold vibration canbe maintained. α and β are natural angular frequencies of thetwo second-order systems, which are known to be almost constant within an utterance as well as across utterances for a particular speaker. It has been shown that α = 3 rad/s and β = 20
rad/s can be used as default values [1].

Given the state sequence s = {sk }K
k=1 , the above HMM
generates the up [k] and ua [k] pair. From (2) and (4), up [k]and ua [k] are then fed through different critically damped ﬁlters, Gp [k] and Ga [k], to generate the phrase and accent components, xp [k] and xa [k]:xp [k] = up [k] ∗ Gp [k],xa [k] = ua [k] ∗ Ga [k],2.2. Probabilistic speech F0 contour model
Here, we brieﬂy review our probabilistic pitch contour model
based on the discrete-time version of the Fujisaki model [5, 6,
7].
In the original Fujisaki model, the phrase commands andaccent commands are assumed to consist of Dirac deltas andrectangular pulses, respectively. In addition, they are not allowed to overlap each other. To incorporate these requirements,
we ﬁnd it convenient to model the up [k] and ua [k] pair, i.e.,o[k] = (up [k], ua [k])T , using a hidden Markov model (HMM).
In [5, 6, 7], we have assumed that {o[k]}K
k=1 is a sequence of
outputs generated from an HMM with the speciﬁc topology illustrated in Fig. 2. The output distribution of each state is assumed to be a Gaussian distribution
o[k] ∼N (o[k]; csk , Υsk ) ,

(7)
(8)

where ∗ denotes convolution over k. An F0 contour is thengiven byx[k] = xp [k] + xa [k] + xb ,(9)where xb denotes the baseline value.
For real speech F0 contours, we must take account of
the uncertainty in the observed F0 data, since observed data
should not always be considered reliable. For example, observed F0 values in unvoiced regions must not be trusted. To
incorporate the degree of uncertainty of F0 observations, weconsider modeling an observed F0 contour y[k] as a superposition of the `“ideal” F´0 contour x[k] and a noise component
xn [k] ∼ N 0, υn2 [k] , where υn2 [k] represents the degree of
uncertainty of the F0 observation at time k, which is assumed
to be given.
Overall, an observed F0 contour y[k] is described as(6)

where sk indicates the state variable. Namely, the mean vector—[k] = (µp [k], µa [k])T = csk and covariance matrix Σ[k] =Υsk are considered to evolve in time as a result of the statetransition s1 , . . . , sK . The deﬁnition of the above HMM can besummarized as follows:y[k] = x[k] + xn [k].

1018

(10)

By marginalizing xn [k] out, we obtain the probability densityKfunction of y = {y[k]}K
k=1 , given o = {o[k]}k=1 , as
P (y|o) =KYN (y[k]; x[k], υn2 [k]),k=1

x[k] = Gp [k] ∗ up [k] + Ga [k] ∗ ua [k] + ub .

(11)

Recall from (6) that given a state sequence s={sk }Ko is generated according to P (o|s, „) =k=1 ,QK
k=1 N (o[k]; csk [k], Υsk ) where „ denotes the mean vectors and covariance matrices of the state emission distributions.
P (s) is given by the product of the state transition probabilities:
Q
P (s) = ϕs1 K
k=2 ϕsk ,sk−1 .

3. Proposed model
Figure 3: Proposed HMM topology for phrase-accent commandsequence modeling based on the vocabulary model of pitch pattern templates.

3.1. Vocabulary model of phrase-accent command sequence
A key idea in the above model is that the sequence of the phrase
and accent command pair is modeled as an HMM. In this paper,
we are concerned with designing an HMM topology (a statetransition network). The previous HMM topology shown in 2.2
allows for any state transitions unless they break the Fujisaki
model’s constraint. However, most of the state transitions are
not appropriate in a linguistic sense. If we can effectively constrain the solution space of the phrase and accent command sequence through an HMM topology design, we may expect toobtain a better solution to the inverse problem. Here, a questionis how to design a state transition network for our model.
In normal speech, we use many intonation patterns when
we speak. However, the types of intonation patterns are limited. We use a common intonation pattern to utter a large set ofphrases (or sentences). As for Japanese, phrases or sentences
can be categorized according to accent types. Phrases with thesame accent type show the same pitch pattern. For example,
when we utter the phrases “arayuru genjitsu wo” and “ashita
wa rinko da”, the intonation patterns should be almost the same
since the corresponding accent types are the same. Thus, it may
be natural to hypothesize that all phrase and accent commandsequences are drawn from a vocabulary consisting of relatively
small number of intonation pattern templates. Namely, we assume that we have a dictionary consisting of a ﬁnite number ofleft-to-right HMM templates and a sequence of the phrase andaccent command pair is generated according to a concatenationof those templates. This can be modeled by an HMM topology
described in Fig. 3. In order to parameterize the duration of theself-transition (except for state p1 , p2 , . . .), as with [6, 7], we
split each state into a certain number of substates such that they
all have exactly the same emission densities. Fig. 4 shows anexample of the splitting of state a1,1 .
The proposed method consists of two stages: training andrecognition. In the training stage, the mean vectors of the state
emission distributions and the state transition probabilities arelearned from a training data set. This corresponds to learning
the intonation pattern templates. In the recognition stage, we
search for the optimal state sequence with ﬁxed templates andstate transition probabilities, given a test data. Both stages are
based upon the same optimization algorithm, which will be decribed in the next section. The only difference is that the HMM
parameters (the mean parameters of the state emission distributions and the state transition probabilities) are ﬁxed during therecognition stage.

Figure 4: The splitting of state a1,1 into 4 substates a1,1,0 ,
a1,1,1 , a1,1,2 , and a1,1,3 . ϕa1,1,0 ,a1,1,1 corresponds to the probability of staying at state a1,1 with 4 consecutive times. states r
are treated similarly with states a.

3.2. Parameter Optimization Algorithm
In this section, we describe an iterative algorithm that searches
for the maximum a posteriori estimates of o and „ by locally maximizing P (o, „|y) given y using the generalized
Expectation-Maximization (EM) algorithm. We treat s asa latent variable and consider marginalizing P (o, „, s|y) ∝P (y|o)P (o|s, „)P (s) with respect to s to obtain the objective P (o, „|y). The auxiliary function (as known as the “Qfunction”) can be written as
X
Q(o, „, o′ , „ ′ ) =P (s|y, o′ , „ ′ ) log P (o, „, s|y)sc

= log P (y|o) +X

P (s|y, o′ , „ ′ ) log P (o|s, „)P (s),scwhere = denotes equality up to constant terms. An iterativealgorithm that consists of computing P (s|y, o′ , „ ′ ) (via the
Forward-Backward algorithm), increasing Q(o, „, o′ , „ ′ ) withrespect to o and „, and then substituting o and „ into o′ and„ ′ locally maximizes the posterior P (o, „|y). Here, care must
be taken that increasing Q(o, „, o′ , „ ′ ) with respect to o mustbe performed subject to non-negativity. This can be done by
invoking the idea of [8]. By using the Jensen’s inequality weobtain an inequality
−

“

X

X

i∈{p,a,b}l

Gi [k − l]ui [l]
≥−

X

”2
X G2i [k − l]u2i [l]i∈{p,a,b} l1019

λi,k,l,(12)

Table 1: Accuracy rates (S=0.3s). The left, middle, and rightcolumns show the accuracy rates of the phrase and accent commands, the phrase commands alone, and the accent commands
alone, respectively. The “Init” row shows the accuracy rates
obtained with the θ before the EM iterations (henceforth Narusawa’s method [4]), the “N templates” row shows that of theestimated command sequence after the EM iterations with thenumber N of the templates, and the “Conv.” row shows thoseobtained with our previous method [6, 7].

0.6

magnitude

Template A

Template B

Template A

Phrase commands
Accent commands

0.4

0.2

0

0

50

100

150

200

250

300

350

400

450

350

400

450

time
0.8
Template B

0.6

Template B

0.4

Init [4]T =5
T = 10
T = 15Conv. [6, 7]

All

Phrase commands

Accent commands

65.9%
70.2%
70.0%
67.3%
69.3%

67.4%
70.3%
71.0%
72.5%
63.8%

65.1%
70.0%
69.5%
64.7%
72.1%

0.2
0

0

50

100

150

200

250

300

time

Figure 5: An example of results from Viterbi decoding in therecognition stage. Commands sequences are extracted fromjapanese utterances “watashi wa sore o ryokan ni motte kaetta”(above) and “ippyo no kakusa wa sarani hirogaru daro” (below.) In the ﬁrst command sequence, there are three templates.
First and second templates are the same and are classiﬁed as“template A,” while the third is classiﬁed as “template B.” In thesecond command sequence, there are two templates and bothare classiﬁed as “template B.”

where Gb [k] = δ[k] (Kronecker’s
P P delta), λi,k,l ≥ 0 is an auxiliary variable satisfying i l λi,k,l = 1. We can use this inequality to construct a lower bound function for Q(o, „, o′ , „ ′ ).
The maximization of this lower bound function w.r.t. o (subject
to non-negativity) and λ can be achieved analytically, whichguarantees a certain increase of Q(o, „, o′ , „ ′ ).
After convergence, we search for the optimal state sequences by using the Viterbi algorithm.

S = 0.3 s. The accuracy rate is slightly improved in T =5 and 10. This is probably because the present method is atemplate-based method which made a best use of information
about the levels and timings of accent commands to accurately
estimate at which position each phrase command should appear.
In T = 15, although the accuracy rate of the phrase commands
was improved, unnecessary templates with 4 accent commands
seemed to make the accuracy of accent commands worse.

4. Experimental Evaluation
Accuracy of phrase and accent commands estimation should depend on the size of the vocabulary model (i.e., the number of theintonation pattern templates). Thus, we evaluated the present
method with different settings of the number of the templates interms of estimation accuracy.
In the training stage, For each utterance in the ﬁrst 50 sentences of ATR speech database [9] spoken by a male speaker(MHT), F0 contours were extracted using the method describedin [10], from which the model parameters were estimated. The
constant parameters were ﬁxed respectively at t0 = 8 ms,
α = 3.0 rad/s, β = 20.0 rad/s, υn2 [k] = 1015 for unvoiced
regions and υn2 [k] = 0.22 for voiced regions. ub was set atthe minimum log F0 value in the voiced regions. Numbers of
substates for states a and a were 250.
In the recognition stage, for each F0 contour extracted in
the same way as training stage from the last 53 sentences in the
ATR speech database, the Fujisaki model parameters were estimated from which commands estimation accuracies were calculated. In both stage, the initial values of θ were set at thevalues obtained with the method described in [4]. We tested
the present method with different settings of the number of thetemplates T : 5, 10, and 15. The number of accent commands ineach template was set at 1 to 3. The numbers of templates with
1, 2, 3 or 4 accent commands were (2, 2, 1, 0), (4, 4, 2, 0) and(5, 5, 3, 2) for T = 5, 10, 15, respectively.
Evaluation procedure is same as our previous works [6, 7],
where a pair of commands are matched in the dynamic programming if the difference between the estimated and groundtruth was shorter than S seconds. The magnitudes of the phrase
and accent commands not taken account because the magnitude
estimation was very sensitive to the baseline F0 value, whichwere set differently in the present method and in the manual
annotation.
Tab. 1 shows the result of our quantitative evaluation with

5. Discussion
The proposed method is a template-based method, thus it can beconsidered that the method classiﬁes commands sequences intoﬁnite kind of templates, as shown in Fig. 5. This side effects
suggests possibility to combine our method with clusteringbased approach, namely speaker recognition [11], voice conversion [12], or text to speech synthesis [13]. Such classiﬁcation can worsen accuracies of command estimation, however,as shown in the preceding sections, accuracy rates are at least
equivalent to previous methods.

6. Conclusion
This paper dealt with the problem of estimating the Fujisakimodel parameters (phrase-accent command sequence) from a
raw F0 contour, and proposed a probabilistic generative modelof speech F0 contours incorporating a vocabulary model ofpitch pattern templates. The present method was evaluated interms of phrase and accent commands estimation accuracy withdifferent settings of the number of the templates. The experimental results revealed that the reduction in the solution spaceof the phrase and accent sequence had an effect on improving
the accuracy of command estimation.

7. Acknowledgement
We thank Prof. Keikichi Hirose of the University of Tokyo, whokindly provided us with the manually annotated ground truth
data associated with the ATR speech samples.

1020

8. References[1] H. Fujisaki, In Vocal Physiology: Voice Production, Mechanisms
and Functions, Raven Press, 1988.
[2] H. Fujisaki and K. Hirose, “Analysis of voice fundamental frequency contours for declarative sentences of japanese,” J. Acoust.
Soc. Jpn (E), vol. 5, no. 4, pp. 233–242, 1984.
[3] H. Mixdorf, “A novel approach to the fully automatic extractionof fujisaki model parameters,” in Proc. ICASSP, 2000, vol. 3, pp.
1281–1284.
[4] S. Narusawa, N. Minematsu, K. Hirose, and H. Fujisaki, “A
method for automatic extraction of model parameters from fundamental frequency contours of speech,” in Proc. ICASSP, 2002,
pp. 509–512.
[5] H. Kameoka, J. Le Roux, and Y. Ohishi, “A statistical model ofspeech F0 contours,” in Proc. SAPA, 2010, pp. 43–48.
[6] K. Yoshizato, H. Kameoka, D. Saito, and S. Sagayama, “Statistical approach to fujisaki-model parameter estimation from speech
signals and its quantitative evaluation,” in Proc. Speech Prosody
2012, 2012, pp. 175–178.
[7] K. Yoshizato, H. Kameoka, D. Saito, and S. 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.
[8] H. Kameoka, T. Nakatani, and T. Yoshioka, “Robust speech dereverberation based on non-negativity and sparse nature of speech
spectrograms,” in Proc. ICASSP, 2009, pp. 45–48.
[9] A. Kurematsu, K. Takeda, Y. Sagisaka, S. Katagiri, H. Kuwabara,
and K. Shikano, “ATR japanese speech database as a tool of speech
recognition and synthesis,” Speech Communication, vol. 27, pp.
187–207, 1999.
[10] H. Kameoka, “Statistical speech spectrum model incorporating all-pole vocal tract model and F0 contour generating processmodel,” in Tech. Rep. IEICE, 2010, in Japanese.
[11] F. Soong, A. Rosenberg, L. Rabiner, and B. Juang, “A vector
quantization approach to speaker recognition,” in Proc. ICASSP,
1985, pp. 387–390.
[12] M. Abe, S. Nakamura, K. Shikano, and H. Kuwabara, “Voice
conversion through vector quantization,” in Proc. ICASSP, 1988,
pp. 655–658.
[13] K. Tokuda, T. Yoshimura, T. Masuko, T. Kobayashi, T. Kitamura,“Speech parameter generation algorithms for HMM-based speechsynthesis,” in Proc. ICASSP, 2000, pp. 1315–1318.

1021