Proceedings of CRAC-2001 (Workshop on Consistent & Reliable Acoustic Cues for Sound Analysis), September 2001.

A Predominant-F0 Estimation Method for Real-world Musical Audio Signals:
MAP Estimation for Incorporating Prior Knowledge about F0s and Tone Models
Masataka Goto
“Information and Human Activity”, PRESTO, Japan Science and Technology Corporation (JST). /National Institute of Advanced Industrial Science and Technology (former Electrotechnical Laboratory).
1-1-1 Umezono, Tsukuba, Ibaraki 305-8568, JAPANm.goto@aist.go.jp
Abstract
In this paper we describe a robust method, called PreFEst, forestimating the fundamental frequency (F0) of melody and basslines in monaural audio signals containing sounds of various
instruments. Most previous F0-estimation methods have difficulty dealing with such complex audio signals because they aredesigned for mixtures of only a few sounds. Without assuming the number of sound sources, PreFEst can obtain the most
predominant F0 — corresponding to the melody or bass line
— supported by harmonics within an intentionally-limited frequency range. It estimates the relative dominance of every possible F0 (represented as a probability density function of the F0)
and the shape of harmonic-structure tone models by using the
MAP (Maximum A Posteriori Probability) estimation considering their prior distribution. Experimental results showed that areal-time system implementing this method is robust enough todetect the melody and bass lines in compact-disc recordings.

1. Introduction
What cues can we rely on when we want to enable a system tounderstand music? Most people think that reliable and fundamental information in music is provided by the notes in musicalscores, and most previous systems attempting to understand music have dealt with automatic music transcription, which transforms audio signals into a symbolic note-level representation.
Such a representation, however, is not a good starting point
for understanding musical audio signals. Identifying the names(symbols) of the notes corresponding to the sounds of music isa skill mastered only by trained musicians. Even though many
trained musicians use scores for composing and interpreting music, untrained listeners understand music without mentally representing audio signals as musical scores.
We therefore proposed a research approach, music scene description [1], intended to obtain more fundamental and reliable
descriptions of musical audio signals. We have been working on
obtaining intuitive and basic descriptions — such as the melody
line, bass line [2, 3], and hierarchical beat structure [4, 5, 6]— without identifying musical notes and without segregating
sound sources. This paper focuses on detecting melody and basslines in compact-disc recordings. These lines are fundamental tothe perception of Western music and are useful in various practical applications, such as automatic music indexing for information retrieval.
This paper describes a predominant-F0 (fundamental frequency) estimation method, called PreFEst [2], that can detect
the melody and bass lines in monaural complex mixtures containing simultaneous sounds of various musical instruments. It
has been considered difficult to estimate the F0 in such audio
signals because the number of sound sources in them generally cannot be assumed, because the frequency components of
one sound often overlap the frequency components of simultaneous sounds, and because the F0’s frequency component (thefrequency component corresponding to the F0) is sometimes

very weak or missing (missing fundamental). Most previous F0estimation methods [7, 8, 9, 10, 11], however, assumed that the
input contained just a single-pitch sound with aperiodic noises.
Although several methods for dealing with multiple-pitch mixtures were proposed [12, 13, 14, 15], they assumed the numberof simultaneous sounds and had difficulty dealing with compactdisc recordings.
Advantages of our method are that it does not assume thenumber of sound sources, locally trace frequency components,
or even rely on the existence of the F0’s frequency component.
It basically estimates the F0 of the most predominant harmonic
structure in the input sound mixture, simultaneously taking into
consideration all the possibilities of the F0 and treating the input
mixture as if it contains all possible harmonic structures withdifferent weights (amplitudes). It regards a probability density
function (PDF) of the input frequency components as a weighted
mixture of the harmonic-structure tone models (represented by
PDFs) of all possible F0s and then finds the F0 of the maximumweight model corresponding to the most predominant harmonic
structure.
The method has recently been made more adaptive and flexible by the following three extensions [3]: introducing multiple
types of harmonic-structure tone models, estimating the shape of
tone models, and introducing a prior distribution of the model
shapes and F0 estimates. For example, when prior knowledge
about very rough F0 estimates of the melody and bass lines isavailable, it can be used in the estimation process. These extensions were made possible by the MAP (Maximum A Posteriori
Probability) estimation executed by using the EM (ExpectationMaximization) algorithm.
The following sections first describe the extended PreFEst
in detail and then present experimental results showing that areal-time system based on the PreFEst can detect the melody andbass lines and that the use of a prior distribution of F0 estimatesenables the F0 to be determined more accurately.

2. Predominant-F0 Estimation Method:
PreFEstPreFEst consists of three components, PreFEst-front-end for
frequency analysis, PreFEst-core estimating the most predominant F0, and PreFEst-back-end considering temporal continuityof the F0. Since the melody line tends to have the most predominant harmonic structure in middle- and high-frequency regions
and the bass line tends to have the most predominant harmonic
structure in a low-frequency region, we can estimate the F0s of
the melody and bass lines by applying PreFEst-core with appropriate frequency-range limitation.
2.1. PreFEst-front-end: Forming the Observed Probability
Density Functions
The PreFEst-front-end first uses an STFT-based multirate filter bank (Figure 1) in order to obtain adequate time and fre-

Audio signals

FFT
FFT
FFT
FFT

16 kHz

Deci

8 kHz

2

Deci

4 kHz

Deci

Deci : Decimator

2 kHz

Deci

(LPF (0.45 fs) + 1/2 down-sampler)

3.6-7.2
1.8-3.6

FFTkHz

0.45-0.9

c (1|F,2)c(t)(2|F,2)

Tone Model (m=1)
p(x | F,1,µ(t)(F,1))c(t)(3|F,2)c(t)(1|F,1)c(t)(4|F,2)c(t)(2|F,1)

0-0.45kHz

Figure 1: Structure of the multirate filter bank: The cut-off frequency of the anti-aliasing filter (FIR LPF) in each decimator is0.45 fs , where fs is the sampling rate at that branch. The input
signal is digitized at 16 bit / 16 kHz and is finally down-sampled
to 1 kHz.
BPF for detecting bass line

(t)

kHz

0.9-1.8
kHz1kHzp(x,1 | F,2,µ(t)(F,2))

Tone Model (m=2)
p(x | F,2,µ(t)(F,2))kHzc(t)(3|F,1)c(t)(4|F,1)

F

fundamental
frequency F+1200

F+1902
F+2400x [cent]x [cent]

Figure 3: Model parameters of multiple adaptive tone models.

BPF for detecting melody line1p(x, h|F, m, µ(t) (F, m))= c(t) (h|F, m) G(x; F + 1200 log2 h, Wi ),0 cent16.35 Hz

2400 cent65.41 Hz4800 cent261.6 Hz

7200 cent1047 Hz

9600 cent
4186 Hz

Figure 2: Frequency responses of bandpass filters (BPFs).

quency resolution under the constraint of real-time operation. It
then extracts frequency components by using an instantaneousfrequency-related measure [1, 2] and obtains two sets of the input bandpass-filtered frequency components (Figure 2), one forthe melody line and the other for the bass line. To use statisticalmethods, we represent each of the bandpass-filtered frequency
components as a probability density function (PDF), called an
observed PDF, p(t)
Ψ (x), where t is the time measured in units offrame-shifts (10 msec), and x is the log-scale frequency denotedin units of cents (a musical-interval measurement). Frequency
fHz in hertz is converted to frequency fcent in cents as follows:fHz
fcent = 1200 log2.
(1)
3
440 × 2 12 −5
2.2. PreFEst-core: Estimating the F0’s Probability Density
Function
For each set of filtered frequency components represented as anobserved PDF p(t)
Ψ (x), the PreFEst-core forms a probability density function of the F0, called the F0’s PDF, p(t)F 0 (F ), where F
is the log-scale frequency in cents. We consider each observed
PDF to have been generated from a weighted-mixture model oftone models of all the possible F0s; a tone model is the PDF corresponding to a typical harmonic structure and indicates where
the harmonics of the F0 tend to occur. Because the weights oftone models represent the relative dominance of every possible harmonic structure, we can regard those weights as the F0’s
PDF: the more dominant a tone model in the mixture, the higher
the probability of the F0 of its model.
2.2.1. Weighted-mixture model of adaptive tone models
To deal with diversity of harmonic structure, the PreFEst-core
can use multiple types of harmonic-structure tone models. The
PDF of the m-th tone model for each F0 F is denoted byp(x|F, m, µ(t) (F, m)) (Figure 3) where the model parameterµ(t) (F, m) represents the shape of the tone model. The number of tone models is Mi (1 ≤ m ≤ Mi ) where i denotes the
melody line (i = m) or the bass line (i = b). Each tone model isdefined byp(x|F, m, µ(t) (F, m)) =

Hi

h=1

p(x, h|F, m, µ(t) (F, m)),(2)µ (F, m) = {c (h|F, m) | h = 1, ..., Hi },
(t)

(t)

(3)
(4)

(x−x0 )2
−
2σ 2

1
e
,
(5)
2πσ 2
where Hi is the number of harmonics considered, W2i is the variance of the Gaussian distribution G(x; x0 , σ), and c(t) (h|F, m)determines the relative amplitude of the h-th harmonic component (the shape of tone model) and satisfies
G(x; x0 , σ) = √

Hi

(t)

c (h|F, m) = 1.

(6)

h=1

We then consider the observed PDF p(t)
Ψ (x) to have been
generated from the following model p(x|θ (t) ), which is aweighted mixture of all possible tone models p(x|F, m, µ(t) (F,m)):
p(x|θ(t) ) = Fhi Mi
Fliw(t) (F, m) p(x|F, m, µ(t) (F, m)) dF,

m=1

θ(t) = {w(t) , µ(t) },

(7)
(8)

w(t) = {w(t) (F, m) | Fli ≤ F ≤ Fhi , m = 1, ..., Mi },(9)µ(t) = {µ(t) (F, m) | Fli ≤ F ≤ Fhi , m = 1, ..., Mi }, (10)where Fli and Fhi denote the lower and upper limits of the possible (allowable) F0 range and w (t) (F, m) is the weight of a tonemodel p(x|F, m, µ(t) (F, m)) that satisfies

 Fhi Mi
Fliw(t) (F, m) dF = 1.

(11)

m=1

Because we cannot know a priori the number of sound sources,
it is important that we simultaneously take into consideration all
the possibilities of the F0 as expressed in Equation (7). If we canestimate the model parameter θ(t) such that the observed PDF(t)p(t)Ψ (x) is likely to have been generated from the model p(x|θ ),(t)(t)
w (F, m) can be interpreted as the F0’s PDF pF 0 (F ):p(t)F 0 (F ) =Miw(t) (F, m) (Fli ≤ F ≤ Fhi ).

(12)

m=1

2.2.2. Introducing a prior distribution
To use prior knowledge about F0 estimates and the tone-model
shapes, we define a prior distribution p0i (θ(t) ) of θ(t) as follows:
p0i (θ(t) ) = p0i (w(t) ) p0i (µ(t) ),
(13)
(t)
(t)
(t)
1
e−βwi Dw (w0i ;w ) ,(14)
p0i (w(t) ) =Zw Fh Mi

1 − Fl ii
p0i (µ ) =
e
Zµ
(t)

m=1

(t)

(t)

βµi (F,m) Dµ (µ0i (F,m);µ(t) (F,m)) dF.
(15)
Here p0i (w(t) ) and p0i (µ(t) ) are unimodal distributions: p0i (w(t) )(t)takes its maximum value at w0i(F, m) and p0i (µ(t) ) takes its(t)(t)(F, m) and µ(t)maximum value at µ0i (F, m), where w0i
0i (F, m)(t)(c0i (h|F, m)) are the most probable parameters. Zw and Zµ(t)(t)
are the normalization factors, and βwi
and βµi(F, m) are theparameters determining how much emphasis is put on the maximum value. The prior distribution is not informative (i.e., it is(t)(t)
uniform) when βwi
and βµi(F, m) are 0, corresponding to thecase when no prior knowledge is available. In Equations (14)(t)(t)
and (15), Dw (w0i
; w(t) ) and Dµ (µ(t)
0i (F, m); µ (F, m)) are thefollowing Kullback-Leibler information:(t); w(t) ) =Dw (w0i Fhi Mi
Fliw(t) (F, m)(t)dF,(16)w0i(F, m) log 0i(t)w (F, m)

m=1

In the E-step, Q(θ(t) |θ (t) ) is expressed as
Q(θ (t) |θ(t) ) =−∞m=1 h=1∂∂w (t)

 ∞ 
Hi−∞ h=1

(t)p(t))Ψ (x) p(F, m, h|x, θ(log w(t) (F, m) + log p(x, h|F, m, µ(t) (F, m))) dx
c(t) (h|F, m). (17)
c0i (h|F, m) log 0ic(t) (h|F, m)

2.2.3. MAP estimation using the EM algorithm

(t)
(t)
w0i
(F, m) log
− βwi

∂
∂c(t)

(t)
(F, m)
w0i
w(t) (F, m)+ log c(t) (h|F, m) + log G(x; F + 1200 log2 h, Wi )) dx(t)(F, m) c(t)− βµi0i (h|F, m) log(18)

Because this maximization problem is too difficult to be solved
analytically, we use the Expectation-Maximization (EM) algorithm [16], which is an iterative algorithm successively applying
two steps — the expectation step (E-step) and the maximization
step (M-step) — to compute MAP estimates from incomplete(t)observed data (i.e., from p(t)
Ψ (x)). With respect to θ , each iter(t)
(t)
(t)
ation updates the old estimate θ = {w , µ } to obtain thenew (improved) estimate θ(t) = {w(t) , µ(t) }.
By introducing hidden (unobservable) variables F , m, andh, which respectively describe which F0, which tone model,
and which harmonic component were responsible for generating each observed frequency component at x, we can specify
the two steps as follows:
1. (E-step)Compute the following QMAP (θ(t) |θ(t) ) for the MAP estimation:
QMAP (θ(t) |θ(t) ) = Q(θ(t) |θ(t) ) + log p0i (θ(t) ),(19)
Q(θ(t) |θ (t) ) =

 ∞

−∞

p(t)
Ψ (x)

EF,m,h [log p(x, F, m, h|θ (t) ) | x, θ(t) ] dx,(20)
where Q(θ(t) |θ (t) ) is the conditional expectation of themean log-likelihood for the maximum likelihood estimation.
EF,m,h [a|b] denotes the conditional expectation of a with respect to the hidden variables F , m, and h with the probability
distribution determined by condition b.
2. (M-step)
Maximize QMAP (θ(t) |θ(t) ) as a function of θ(t) in order toobtain the updated (improved) estimate θ(t) :θ (t) = argmax QMAP (θ(t) |θ(t) ).
(21)
θ (t)

= 0, (24)

(t)
p(t)
) (log w(t) (F, m)
Ψ (x) p(F, m, h|x, θ

−∞

∞

(t)
(t)
p(t)
Ψ (x) (log p(x|θ ) + log p0i (θ )) dx.

− λw (w(t) (F, m) − Mi (Fh1i −Fli ) ) ∞The problem to be solved is to estimate the model parameter
θ(t) , taking into account the prior distribution p0i (θ(t) ), when we
observe p(t)
Ψ (x). The MAP (Maximum A Posteriori Probability)
estimatorof θ (t) is obtained by maximizing
−∞

p(t)
Ψ (x)

= log(w (t) (F, m) p(x, h|F, m, µ(t) (F, m))).
(23)
Regarding the M-step, Equation (21) is a conditional problem ofvariation, where the conditions are given by Equations (6) and(11). This problem can be solved by using the following EulerLagrange differential equations with Lagrange multipliers λw
and λµ :

Hi

(t)
h=1

Fli

p(F, m, h|x, θ(t) ) log p(x, F, m, h|θ (t) )dF dx, (22)where the complete-data log-likelihood is given bylog p(x, F, m, h|θ(t) )(t)Dµ (µ(t)0i (F, m); µ (F, m))= ∞  Fhi 
Mi Hi


c(t)
0i (h|F, m)c(t) (h|F, m)− λµ (c(t) (h|F, m) − H1i )= 0. (25)

From these equations we getw(t) (F, m) =∞−∞(t) (t)(t)p(t)) dx + βwi
w0i (F, m)
Ψ (x) p(F, m|x, θ(t)
1 + βwi,(26)

c(t) (h|F, m) =

∞

−∞

(t)
(t)
p(t)
) dx + βµi(F, m)c(t)Ψ (x)p(F, m, h|x, θ
0i (h|F, m)∞−∞(t) ) dx + β (t) (F, m)p(t)µiΨ (x)p(F, m|x, θ.

(27)
According to the Bayes’ theorem, p(F, m, h|x, θ(t) ) is given byw(t) (F, m) p(x, h|F, m, µ(t) (F, m))p(F, m, h|x, θ(t) ) =.
p(x|θ (t) )(28)
Finally we obtain the following new parameter estimates:
w(t) (F, m) =(t)(t) (t)wML(F, m) + βwi
w0i (F, m)(t)
1 + βwi,(29)

c(t) (h|F, m) =
(t)
(t)
(t)
wML
(F, m) c(t)ML (h|F, m) + βµi (F, m)c0i (h|F, m)
(t)
(t)
wML
(F, m) + βµi(F, m), (30)(t)(F, m) and c(t)where wML
ML (h|F, m) are, when the noninforma(t)(t)
tive prior distribution (βwi
= 0 and βµi(F, m) = 0) is given, thefollowing maximum
likelihood
estimates:

∞

(t)
wML
(F, m) =

−∞

p(t)
Ψ (x)

w(t) (F, m) p(x|F, m, µ(t) (F, m)) Fhi MiFli

ν=1w(t) (η, ν) p(x|η, ν, µ(t) (F, ν)) dηdx, (31)c(t)ML (h|F, m) =

1
(t)
wML
(F, m)

 ∞
−∞

p(t)
Ψ (x)

w(t) (F, m) p(x, h|F, m, µ(t) (F, m)) Fhi MiFli

ν=1w(t) (η, ν) p(x|η, ν, µ(t) (F, ν)) dηdx.(32)After the above iterative computation, the F0’s PDF p(t)F 0 (F )estimated by considering the prior distribution can be obtained
from w (t) (F, m) according to Equation (12). We can also obtainthe tone-model shape c(t) (h|F, m), which is the relative amplitude of each harmonic component of all the types of tone modelsp(x|F, m, µ(t) (F, m)).
2.3. PreFEst-back-end:
Sequential F0 Tracking by
Multiple-Agent Architecture
A simple way to identify most predominant F0 is to find the frequency that maximizes the F0’s PDF. This result is not stable,
however, because peaks corresponding to the F0s of simultaneous sounds sometimes compete in the F0’s PDF for a moment
and are transiently selected, one after another, as the maximum.
We therefore consider the global temporal continuity of theF0 by using a multiple-agent architecture [1, 2] in which agents
track different temporal trajectories of the F0. The final F0 output is determined on the basis of the most dominant and stable
F0 trajectory.

3. Experimental Results
The PreFEst has been implemented in a real-time system thattakes a musical audio signal as input and outputs the detected
melody and bass lines in several forms, such as audio signals
and computer graphics [1, 2]. The current implementation uses
two adaptive tone models with the following parameter values:
Fhm = 8400 cent, Flm = 3600 cent, Mm = 2, Hm = 16, Wm = 17
cent; Fhb = 4800 cent, Flb = 1000 cent, Mb = 2, Hb = 6, and Wb= 17 cent. For the prior distribution of the shape of tone models,
we use c(t)0i (h|F, m) = αi,m gm,h G(h; 1, Ui ), where m is 1 or
2, αi,m is a normalization factor, gm,h is 2/3 (when m = 2 andh is even) or 1 (otherwise), Um = 5.5, and Ub = 2.7.
The system was tested on excerpts from a total of 10 songsin the popular, jazz, and orchestral genres. The input monaural
audio signals — each containing a single-tone melody and thesounds of several instruments — were sampled from compact
discs. We evaluated the detection rates by comparing the estimated F0s with the correct F0s hand-labeled using the F0 editor
program we previously developed [2]. In our experiment thesystem correctly detected the melody and bass lines for most ofeach audio sample: the average detection rate was 88.4% for the
melody line and 79.9% for the bass line.
We also tested the system by providing prior knowledge
about rough F0 estimates of the melody line. For the prior F0(t)(F, m) = G(F ; F0i(t) , 100 cent)/Mi ,distribution, we used w0i(t)where F0i is the F0 estimate given by playing a MIDI keyboard
while listening to each excerpt of the songs. Comparing the results obtained with and without using the prior F0 distribution
showed that the use of prior knowledge improved the average
detection rate for the melody line (from 88.4% to 91.2%). In
particular, the detection rate for an orchestral song was greatly
improved (11.9% improvement).

4. Conclusion
We have described a method, called PreFEst, that estimates themost predominant F0 in a monaural complex sound mixture
without assuming the number of sound sources. Its MAP estimation executed by using the EM algorithm makes it possible toestimate the F0’s PDF and the shape of tone models while considering their prior distribution. Experimental results showedthat a system implementing this method is robust enough to estimate the F0s of the melody and bass lines in compact-discrecordings in real time.
Although the PreFEst has great potential, we have not fully
exploited it. In the future, for example, a lot of tone models
could be prepared by analyzing various kinds of harmonic structure appearing in music, and multiple peaks in the F0’s PDF,
each corresponding to a different sound source, could be tracked
simultaneously by using a sound source discrimination method.
The PreFEst can also be applied to non-music audio signals. In
fact, Masuda-Katsuse [17, 18] has extended it and shown that
it is effective for speech recognition in realistic noisy environments.
Acknowledgments: I thank Shotaro Akaho and Hideki Asohfor their valuable discussions.

5. References[1] M. Goto, “A real-time music scene description system: Detecting melody and bass lines in audio signals,” Working Notes of the
IJCAI-99 Workshop on Computational Auditory Scene Analysis,
pp. 31–40, 1999.
[2] M. Goto, “A robust predominant-F0 estimation method for realtime detection of melody and bass lines in CD recordings,” Proc.
ICASSP 2000, pp. II–757–760, 2000.
[3] M. Goto, “A predominant-F0 estimation method for CD recordings: MAP estimation using EM algorithm for adaptive tone models,” Proc. ICASSP 2001, 2001.
[4] M. Goto and Y. Muraoka, “A beat tracking system for acoustic
signals of music,” Proc. the Second ACM Intl. Conf. on Multimedia, pp. 365–372, 1994.
[5] M. Goto and Y. Muraoka, “Real-time beat tracking for drumlessaudio signals: Chord change detection for musical decisions,”
Speech Communication, vol. 27, nos. 3–4, pp. 311–335, 1999.
[6] M. Goto, “An audio-based real-time beat tracking system for music with or without drum-sounds,” Journal of New Music Research, vol. 30, 2001. (in press)[7] L. R. Rabiner et al., “A comparative performance study of several
pitch detection algorithms,” IEEE Trans. on ASSP, vol. ASSP-24,
no. 5, pp. 399–418, 1976.
[8] A. Nehorai and B. Porat, “Adaptive comb filtering for harmonic
signal enhancement,” IEEE Trans. on ASSP, vol. ASSP-34, no. 5,
pp. 1124–1138, 1986.
[9] F. J. Charpentier, “Pitch detection using the short-term phase spectrum,” Proc. ICASSP 86, pp. 113–116, 1986.
[10] T. Abe et al., “Robust pitch estimation with harmonics enhancement in noisy environments based on instantaneous frequency,”
Proc. ICSLP 96, pp. 1277–1280, 1996.
[11] H. Kawahara et al., “Fixed point analysis of frequency to instantaneous frequency mapping for accurate estimation of F0 and periodicity,” Proc. Eurospeech 99, pp. 2781–2784, 1999.
[12] H. Katayose and S. Inokuchi, “The kansei music system,” Computer Music Journal, vol. 13, no. 4, pp. 72–77, 1989.
[13] G. J. Brown and M. Cooke, “Perceptual grouping of musicalsounds: A computational model,” Journal of New Music Research, vol. 23, pp. 107–132, 1994.
[14] K. Kashino and H. Murase, “Music recognition using note transition context,” Proc. ICASSP 98, pp. 3593–3596, 1998.
[15] A. P. Klapuri, “Multipitch estimation and sound separation by thespectral smoothness principle,” Proc. ICASSP 2001, 2001.
[16] A. P. Dempster et al., “Maximum likelihood from incomplete data
via the EM algorithm,” J. Roy. Stat. Soc. B, vol. 39, no. 1, pp. 1–
38, 1977.
[17] I. Masuda-Katsuse, “A new method for speech recognition in thepresence of non-stationary, unpredictable and high-level noise,”Proc. Eurospeech 2001, pp. 1119–1122, 2001.
[18] I. Masuda-Katsuse and Y. Sugano, “Speech estimation biased byphonemic expectation in the presence of non-stationary and unpredictable noise,” Proc. CRAC workshop, 2001.