I-DIVERGENCE-BASED DEREVERBERATION METHODWITH AUXILIARY FUNCTION APPROACH
Naoki Yasuraoka,† Hirokazu Kameoka,‡ Takuya Yoshioka,‡ Hiroshi G. Okuno††‡

Graduate School of Informatics
Kyoto University
Yoshida-honmachi, Sakyo-ku
Kyoto 606-8501, JapanNTT Communication Science Laboratories
NTT Corporation
3-1 Morinosato Wakamiya, Atsugi,
Kanagawa, 243-0198, Japan

ABSTRACTThis paper presents a dereverberation method based on I-divergence
minimization, which is particularly suitable for music signals. Existing dereverberation methods, including one designed for music,
sometimes distort instrument sounds and make staccato-like tones.
The problems with the Itakura-Saito-divergence-based formulationof the existing methods are 1) their tendency to excessive suppression of direct sound and 2) the difﬁculty of incorporating and optimizing sophisticated source models suitable for music signals. The
proposed I-divergence-based method can mitigate these problems.
Employing the I-divergence measure enables us to avoid the direct sound suppression problem and to use powerful music spectrum models without complicating its optimization. We develop
a convergence-guaranteed parameter estimation algorithm based on
the auxiliary function approach. Experimental results reveal the effectiveness of the proposed dereverberation method.
Index Terms— dereverberation, musical audio signal, Idivergence, auxiliary function method
1. INTRODUCTION
Audio effects such as reverberation, delay, and phase shift are important factors for enriching the perceptual quality of music. On
the other hand, they are very harmful for automatic analysis (e.g.
melody extraction [1] and chord detection [2]) of musical audio signals. Therefore, technologies for removing and controlling such audio effects will be vital for development of emerging music processing applications including automatic music analyzers [3] and active
music listening systems [4].
Dereverberation for music signals, namely music dereverberation, is a much more challenging task than speech dereverberation.
This is because musical tones tend to last long unlike moras in natural voices, which makes it difﬁcult to distinguish direct sound fromreverberation. Indeed, existing speech dereverberation methods suchas [5] are not necessarily effective for music dereverberation. Recently, [6] proposed a music dereverberation method that used aharmonicity-based source model by taking advantage of strong harmonic structures inherent in music signals. This method was shown
to be effective in removing reverberation from music signals.
However, we have recently experienced a serious problem evenwith the method of [6]: dereverberated signals with this method, especially for struck string instruments, tend to be distorted and soundlike staccato tones.
We ﬁnd the main cause of this problem to be that this method(as well as other speech dereverberation methods) is based on the

978-1-4577-0539-7/11/$26.00 ©2011 IEEE369

minimization of the Itakura-Saito (IS) divergence. The main drawback of employing the IS-divergence is twofold. Firstly, with the
IS-divergence, we sometimes regard sustained parts of musical tones
as reverberation and suppress the sustain mistakenly, which resultsin the staccato-like sounds. Secondly, because of technical inconveniences in the analysis of the IS-divergence, the parameters of anassumed source model are difﬁcult to estimate. To avoid the latter,[6] settled for using different criteria for dereverberation ﬁlter optimization and source model parameter estimation. But, because ofthis settlement, [6] had to use a non-convergent iterative algorithm.
This paper presents a new dereverberation method based on the
I-divergence criterion allowing for the joint estimation of dereverberation ﬁlter coefﬁcients and source model parameters. In audio processing, the I-divergence is often used as a goodness of ﬁt measure
between an observed spectrum and a source model. The I-divergencemeasure allows us to use powerful music spectrum models such as
harmonic Gaussian mixture model [7, 8] and nonnegative matrix factorization (NMF) model [9, 10]. To realize the I-divergence-based
dereverberation, a new convergence-guaranteed optimization algorithm is proposed based on the auxiliary function concept (e.g., [10]).
The remainder of this paper is organized as follows: Section
2 describes the dereverberation problem and qualitatively compares
existing and proposed reverberation methods. Section 3 introduces
the proposed dereverberation method with I-divergence. Section 4
presents the evaluation results, and Section 5 concludes the paper.
2. COMPARISON OF EXISTING AND PROPOSEDDEREVERBERATION METHODSThis section ﬁrst deﬁnes the dereverberation problem. Then a new
I-divergence-based dereverberation approach is introduced and compared with the existing IS-divergence-based one.
2.1. Problem Statement
Let sn,l denote an anechoic signal of speech or music in the shorttime Fourier transform (STFT) domain, where, n and l correspondto the time frame and frequency bin indices, respectively. We refer tosn,l as a source signal. The source signal is assumed to be unobservable, and only its reverberated version, denoted by xn,l , is available.
Dereverberation is the process of estimating the source signal, sn,l ,by using the reverberant signal, xn,l , without precise knowledge ofthe reverberation properties. We assume that the observed signal
is monaural and dereverberation is performed in a batch processing
manner.
We adopt the inverse ﬁltering approach to achieve dereverberation. As in [5], we assume that the reverberant signal is generated

ICASSP 2011

from the source signal sn,l by an auto regressive (AR) system as
xn,l =

Kg d,l xn−d,l + sn,l .

Table 1. Algorithm for parameter estimation
(1)

d=k

k and K are the initial delay time and the ﬁlter length, respectively.
Then the source signal can be recovered as xn,l − K
d=k g d,l xn−d,l .
Thus, our goal is to estimate g d,l . In this paper we refer to g d,l as areverberation ﬁlter coefﬁcient.
In order to estimate g d,l , we introduce a source model characterized by a set of parameters. The source model is used to measure
the degree to which a dereverberated signal is good. Since the source
model parameters are unavailable in advance, they are also estimatedin the course of the reverberation ﬁlter estimation.
To put this idea into shape, we have to set up an optimization
problem, which is speciﬁed by the following two components:
1. source model reﬂecting the source characteristics such as harmonicity, which is described in Section 3.2, and
2. cost function between the dereverberated signal, rn,l =
xn,l − K
d=k g d,l xn−d,l , and this source model.
2.2. IS-divergence-based existing method
Yoshioka et al. proposed a dereverberation method [5], which minimizes a cost function deﬁned as follows:N−1 L−1|rn,l |2log pn,l +,(2)Q(IS) =
pn,l
n=0
l=0

where pn,l is the source model that deﬁnes the power spectral envelope of the source. From the viewpoint of the estimation of pn,l , thiscost function consists of the Itakura-Saito (IS) divergence between|rn,l |2 and pn,l . Although it is not the IS-divergence in terms of theestimation of g d,l , we simply call it the IS-divergence. Note that thecost function is quadratic form with respect to the ﬁlter coefﬁcient
g d,l , and thus g d,l can be estimated easily.
However, minimizing Eq.(2) with respect to g d,l sometimes
leads to excessive suppression of direct sound. This can be accounted for as follows. We note that rn,l is the prediction error of
the AR system representing reverberation (see Eq.(1)), and Eq.(2) isa weighted sum of the errors. Hence, this cost function tends to confuse the sustain part of direct sound as reverberation because the sustain part is almost predictable from its preceding part. Furthermore,sophisticated source models that were previously proposed for music
spectrum modeling are difﬁcult to optimize with the IS-divergence.
2.3. Proposed I-divergence-based method
In this paper, we consider the I-divergence between the magnitudeof dereverberated signal |rn,l | and the source model pn,l as:N−1 L−1|rn,l |
− (|r n,l | − pn,l ) . (3)
Q(I) =|rn,l | log
pn,l
n=0
l=0

Note that pn,l represents in turn the magnitude spectrum of thesource. Advantages of the proposed approach are twofold. Firstly,with this cost function g d,l is estimated in a way that directly minimize the difference between the dereverberated signal and the source
estimate, which avoids sound distortion. Secondly, various magnitude spectrum modeling techniques for music signals based on the
I-divergence (e.g., Harmonic GMM [7, 8] designed for F0 estimation) can be incorporated into the present dereverberation system.
It is important to note here that the cost function given by Eq.(3)is no longer a quadratic form with respect to the reverberation ﬁlter

370

1. Initialize source magnitude model parameters assuming |rn,l |to be |xn,l |, and initialize g d,l by the IS-divergence based dereverberation given by Eq.(2) with pn,l = |xn,l |2 .
2. Update source magnitude model parameters using current g d,l .
3. Update g d,l using a current source magnitude model estimate:
3.a Update auxiliary variables {tn,l , bn,l , ωn,l } by Eq.(13).
3.b Update g d,l by Eq.(14).
3.c Iterate from 3.a for several times.
4. Iterate 2. and 3. until convergence.
coefﬁcients g d,l . We are thus concerned with deriving an algorithm
for estimating g d,l , as will be discussed in the following section.
3. DEREVERBERATION WITH I-DIVERGENCE
This section describes the proposed parameter estimation algorithm
using the auxiliary function method. The principle of the auxiliary
function method is to iteratively minimize an upper bound of a costfunction, which is easier to minimize than the original cost function.
The upper bound is called an auxiliary function, and it is generally
characterized by both the parameters of interest and additional parameters called auxiliary parameters. See [10] for more details.
3.1. Parameter estimation algorithm
Now our objective is to introduce an iterative algorithm for estimating g d,l efﬁciently on the basis of I-divergence. The key to the development of the proposed algorithm is its use of the auxiliary function
method to deduce from Eq.(3) a quadratic form with respect to g d,l ,which can be optimized analytically. The auxiliary function method
is used twice before obtaining the quadratic form. The overall parameter estimation algorithm is summarized in Table 1.
First we introduce one inequality to eliminate nonlinearity of alogarithmic function in the ﬁrst term in Eq.(3) as follows:1(|rn,l | − tn,l ),(4)
log |rn,l | ≤ log tn,l +
tn,lwhere tn,l > 0 is an auxiliary parameter. The right-hand side of theinequality is a linear function that is tangent to log |rn,l | at |rn,l | =tn,l , and thus can be considered as an auxiliary function of log |rn,l |.
By using this, we obtain an auxiliary function Q+ for Q(I) as:N−1 L−1 |rn,l |2 
tn,l
+
Q =
+ log
− 2 |rn,l | + pn,l . (5)
tn,l
pn,l
n=0
l=0

However, ∂g∂d,l Q+ = 0 still cannot be solved analytically becauseof the second term in Eq.(5).
Next, we introduce the following two inequalities, which arechosen depending on log(tn,l /pn,l ) − 2, appearing in Eq.(5).
in case Dn,l ≥ 0 :in case Dn,l < 0 :|rn,l |2bn,l+, bn,l > 0
2bn,l
2
 ∗
|rn,l | ≥ Re ωn,l rn,l , |ωn,l | = 1|rn,l | ≤(6)(7)

where Dn,l = log(tn,l /pn,l ) − 2 and superscript ∗ stands for complex conjugate. bn,l and ωn,l are auxiliary parameters for Eq.(6) andEq.(7), respectively. The right-hand side of the ﬁrst inequality is aconvex quadratic function that is tangent to |rn,l | at |rn,l | = bn,l .
The right-hand side of the second inequality means inner product ofr n,l and an directional vector ωn,l .

By introducing the above two inequalities, we obtain anothernew auxiliary function Q++ :Q++ =N−1 L−1



Q++
n,l ,

(8)

n=0 l=0


⎧ 
1
Dn,l
⎪
⎪
|rn,l |2 + c
+
⎨ 2b
tn,l
n,l
++

2
Qn,l =
⎪⎪⎩ 1 rn,l + tn,l Dn,l ωn,l  + c

tn,l 
2(Dn,l ≥ 0)
, (9)(Dn,l < 0)tn,lThen Eq.(8) can be rewritten as

2
K


N
−1 L−1κn,l − d=k g d,l xn−d,l 


Q++ =
+ c.
ρn,l
n=0

(Ω −jμi,n )2
− l
2σ 2

, I m,n,l = e−(Ωl −υ m )2
2γ 2

,

(18)

where ui,j,n is the intensity of j-th peak of the i-th tone, μi,n isthe fundamental frequency, and σ 2 is the spectral spread of eachharmonic component. Ωl means a scaling function that maps theindex of frequency bins to Hertz. Inharmonic components such as
drum sound and noise ﬂoor, are modeled by an inharmonic GMM
in the second term of Eq.(17), which has M Gaussians with ﬁxed
means υ m and ﬁxed large variance γ 2 .
Parameters to be estimated are {ui,j,n , v m,n , μi,n }. Below are
update rules of these with I and IS divergence, which are also basedon the auxiliary function method.
(IS)• Let us put auxiliary parameters, Ψi,j,n,l and Ψi,j,n,l , asui,j,n H i,j,n,l|rn,l | , and
(19)
pn,l


ui,j,n H i,j,n,l (|rn,l |2− pn,l )(IS)
Ψi,j,n,l = ui,j,n H i,j,n,l 1 +. (20)
pn,l 2(I)This form is obviously quadratic with respect to g d,l . Thus g d,l canbe estimated in the same way as IS-divergence-based estimation.
g d,l is estimated by alternatively minimizing Q++ with respect
to the auxiliary parameters, {tn,l , bn,l , ωn,l }, and reverberation ﬁlter
coefﬁcients g d,l . The update rules for the auxiliary parameters are(13)

The ﬁrst and the second rules are given by minimizing the righthand side of Eqs. (4) and (6), respectively. The third rule is given by
maximizing the right-hand side of Eq.(7). The update rule for g d,l isderived based on Eq.(12) and we obtain[g k,l , ..., g K,l ]T = Rl−1 rl

m=1

j=1

(I)

(12)

l=0

bn,l = |rn,l | , and ωn,l = rn,l / |rn,l | .

i=1

H i,j,n,l = ewhere c represents terms that do not depend on g d,l . For simplicity,
let us set κn,l and ρn,l as(Dn,l ≥ 0)xn,lκn,l =
, and
(10)
tn,l
xn,l + 2 Dn,l ωn,l (Dn,l < 0)⎧Dn,l
⎨ t 1 + 2b
(Dn,l ≥ 0)
1
n,l
n,l
.
(11)
=
1
⎩
ρn,l
(Dn,l < 0)tn,l = |rn,l | ,spectrum of musical tones by using a GMM where the means ofeach Gaussian component appear at the harmonic frequencies of thetone. Assuming I musical tones, each of which has J harmonics,
are present at each time frame, the magnitude or power of source ismodeled as
 J

I
M



ui,j,n H i,j,n,l +v m,n I m,n,l ,(17)pn,l =(14)

where, Rl and rl are called a modiﬁed correlation matrix and a modiﬁed correlation vector, respectively, and are structured as follows:⎤⎡ N −1N−1 xn−k,l ∗
xn−K,l ∗xn−k,l · · ·xn−k,l ⎥
⎢
⎥
⎢
ρn,l
ρn,l
n=0
⎥
⎢ n=0
.
.
..
⎥ , (15).
.
Rl = ⎢.
⎥⎢.
.
⎥⎢ N −1N−1⎦⎣  xn−k,l ∗ xn−K,l ∗xn−K,l · · ·xn−K,l
ρn,l
ρn,l
n=0
 Nn=0
T
−1
N
−1


κn,l ∗
κn,l ∗rl =. (16)
xn−k,l···xn−K,l
ρn,l
ρn,l
n=0
n=0
Note that Rl and rl can be calculated efﬁciently using the Fast
Fourier Transform in the same way as [11].

Ψi,j,n,l =• Update rules of ui,j,n are:⎧ L−1 (I)⎪l=0 Ψi,j,n,l⎪⎪(I-divergence)⎪ L−1⎨l=0 H i,j,n,l.
ui,j,n =(IS)
L−1⎪⎪ 1  Ψi,j,n,l⎪⎪(IS-divergence)
⎩ L
H
l=0

(21)

i,j,n,lUpdate rules of v m,n are similar to those of ui,j,n .
• An update rule of μi,n with I-divergence is:J L−1(I)
j=1
l=0 Ωl jΨi,j,n,l.
μi,n = J L−1
(I)
2
j=1
l=0 j Ψi,j,n,l(22)

μi,n is hard to update analytically with IS-divergence and obtained in numerical calculations such as the Newton method.
3.2.2. Nonnegative-Matrix-Factorization-based model (NMF)
NMF is a dimension reduction technique that has been successfully
applied to musical signal processing such as source separation, transcription [9] and so on. Update rules of NMF with I and IS divergence are described in [12].
4. EXPERIMENTS

3.2. Source modelThis section introduces two kinds of source models that have estimation algorithms with both I and IS divergence. Convergence isguaranteed in case these models are estimated with the same divergence as one used for estimation of reverberation ﬁlter coefﬁcients.

The accuracy of dereverberation and that of source modeling areevaluated by using simulated data to compare proposed and existing methods.
4.1. Experimental Condition

3.2.1. Harmonic GMM with residual (HGMM)
HGMM is a source model that can represent harmonic structure of
instrument sounds. The HGMM describes a magnitude or power

371

The proposed I-divergence-based dereverberation method is compared with the existing IS-divergence-based method using two different source models, the HGMM and the NMF model. A combina-

4.3kHz

Table 2. Experimental results: The left-hand side in each cell shows the LSD
between the true source and the dereverberated signal. The right-hand is the
LSD between the true source and the signal given by the source model. Smaller
values indicate higher quality. The best performing values are shown in boldface.
Derev.-I / HGMM-I(Proposed)
Derev.-I / NMF-I(Proposed)
Derev.-IS / HGMM-I [6]Derev.-IS / HGMM-IS
Derev.-IS / NMF-IS

k=3
k = 10k = 20k=3
k = 10k = 20k=3
k = 10k = 20k=3
k = 10k = 20k=3
k = 10k = 20K = 80
8.72
8.68
8.65
8.70
8.66
8.68
8.74
9.06
8.64
8.98
8.64
8.78
9.82
9.66
9.07
8.98
8.79
8.76
10.06 13.53
9.24 12.69
8.91 13.37
10.14 10.06
9.24
9.80
8.90
9.16

K = 160
8.68
8.71
8.63
8.72
8.63
8.66
8.78
8.96
8.62
9.08
8.62
8.78
9.98
9.82
9.14
9.19
8.82
8.78
10.27 13.61
9.34 12.80
8.98 13.38
10.34 10.21
9.36 10.15
8.94
9.22

K = 240
8.75
8.78
8.68
8.69
8.62
8.65
8.88
9.03
8.64
8.82
8.63
8.78
10.05
9.89
9.16
9.09
8.85
8.82
10.34 13.64
9.32 13.37
9.01 13.41
10.40 10.26
9.22
9.32
8.96
9.10

3.2
2.1
1.1
0.0

(a) Source signal
4.3kHz
3.2
2.1
1.1
0.0

(b) Reverberated signal
4.3kHz
3.2
2.1
1.1
0.0

(c) Dereverberated signal obtained with the proposed method
4.3kHz
3.2
2.1
1.1
0.0

(d) Dereverberated signal obtained with the existing method

Fig. 1. Example of a dereverberated signal
tion of IS-divergence-based dereverberation and I-divergence-based
HGMM, which is our conventional method [6], is also evaluated.
Twelve musical pieces in RWC Music Database, Jazz and Classic [13] synthesized with a MIDI tone generator, were used as thetrue source signals. The signals were truncated up to 20 seconds,
and reverberated by being convolved with two impulse responses.
The reverberation times of these impulse responses were about 1
and 1.5 seconds, respectively. We measured the dereverberation performance using the log spectral distance (LSD), which is deﬁned as:N−1 L−1
20 log |y n,l |  N L,
(23)
LSD(Y, Z) ≡10|z n,l | n=0
l=0

We calculated LSD between 1) the true source, S = {sn,l }, and the
dereverberated signal, R = {r̂n,l }, and 2) the true source and thesignal given by the source model, P = {p̂n,l }, where ˆ means theestimated result. The former, LSD(S, R), stands for dereverberation
accuracy, and the latter, LSD(S, P), stands for source modeling accuracy. Sampling rate was 44.1kHz. STFT is computed using a 1024
point Gaussian window that has a 256-point overlap.
4.2. Results and Discussion
Table 2 shows the evaluation results with several conditions about
initial delay time k and ﬁlter length K. Each left-hand value in eachcell shows LSD(S, R) and each right-hand value shows LSD(S, P).
First, we can see that the accuracy of the I-divergence-based source
model estimation was better than that of the IS-divergence-basedsource model estimation. Furthermore, if we compare the resultsof I-divergence-based dereverberation and the method based on [6],we ﬁnd that this improvement in source model estimation led to the
dereverberation performance improvement only when the reverberation ﬁlter is also estimated based on the I-divergence. This clearly
indicates the advantage of the use of the I-divergence both for source
model estimation and reverberation ﬁlter estimation.
Figure 1 shows the piano spectrograms of the dereverberated
signal obtained with the proposed (c) and existing (d) methods. It
can be found that (c) could successfully cancel the reverberation andavoid distorting the sustain of piano tones.
Another small experiment for real reverberation on CDs showedsimilar property as the above results. The experiment also showedthat the dereverberation performance did not vary much according
to the signal length. The audio ﬁles of both simulated data and real
recordings are available at http://winnie.kuis.kyoto-u.
ac.jp/members/yasuraok/ICASSP2011/.

372

5. CONCLUSION
This paper presented a new dereverberation method using Idivergence, which is often used for musical signal processing.
The proposed method uses the auxiliary function method forconvergence-guaranteed parameter estimation. Our experimental
results revealed that the proposed method could perform music
dereverberation more accurately and robustly than the existing ISdivergence-based one.
Future work includes integrating the proposed method into various types of musical signal processing such as sound source separation and sound analysis and synthesis. We will also extend our
method for multi-channel processing, in particular, stereo input.
6. REFERENCES[1] A. P. Klapuri, “Multiple fundamental frequency estimation based on
harmonicity and spectral smoothness,” IEEE Trans. Speech and Audio
Process., vol. 11, no. 6, pp. 804–816, 2003.
[2] K. Lee and M. Slaney, “Acoustic chord transcription and key extraction
from audio using key-dependent HMMs trained on synthesized audio,”IEEE Trans. Audio, Speech and Lang. Process., vol. 16, no. 2, pp. 291–
301, 2008.
[3] J. A. Moorer, On the segmentation and analysis of continuous musicalsound by digital computer., Ph.D. thesis, Stanford University, 1975.
[4] M. Goto, “Active music listening interfaces based on signal processing,” in Proc. ICASSP, 2007, vol. 4, pp. 1441–1444.
[5] T. Yoshioka, T. Nakatani, and M. Miyoshi, “An integrated method forblind separation and dereverberation of convolutive audio mixtures,” inProc. EUSIPCO, 2008.
[6] N. Yasuraoka, T. Yoshioka, T. Nakatani, A. Nakamura, and H. G.
Okuno, “Music dereverberation using harmonic structure source modeland Wiener ﬁlter,” in Proc. ICASSP, 2010, pp. 53–56.
[7] H. Kameoka, T. Nishimoto, and S. Sagayama, “Extraction of multiple
fundamental frequencies from polyphonic music using harmonic clustering,” in Proc. ICA, 2004, pp. I–59–62.
[8] H. Kameoka, T. Nishimoto, and S. Sagayama, “A multipitch analyzer
based on harmonic temporal structured clustering,” IEEE Trans. Audio,
Speech and Lang. Process., vol. 15, no. 3, pp. 982–994, 2007.
[9] P. Smaragdis and J.C Brown, “Non-negative matrix factorization for
polyphonic music transcription,” in Proc. WASPAA, 2003, pp. 170–180.
[10] D. D. Lee and H. S. Seung, “Algorithms for non-negative matrix factorization,” in Proc. NIPS, 2001, pp. 556–562.
[11] T. Yoshioka, T. Nakatani, and M. Miyoshi, “Fast algorithm for conditional separation and dereverberation,” in Proc. EUSIPCO, 2009.
[12] M. Nakano, H. Kameoka, J. Le Roux, Y. Kitano, N. Ono, and
S. Sagayama, “Convergence-guaranteed multiplicative algorithms fornonnegative matrix factorization with beta-divergence,” in Proc. MLSP,
2010.
[13] M. Goto, H. Hashiguchi, T. Nishimura, and R. Oka, “RWC music
database: Popular, classical, and jazz music databases,” in Proc. ISMIR,
2002, pp. 287–288.